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In this thesis wc study electron transport through magnetic nanocontacts and nanowires 
with ab initio quantum transport calculations. The aim is to gain a thorough understand- 
ing of the interplay between electrical conduction and magnetism in atomic-size conduc- 
tors and how it is affected by different aspects as e.g. the atomic structure and the chem- 
ical composition of the conductor. To this end our ab initio quantum transport program 
ALACANT which combines the non-equilibrium Green's function formalism (NEGF) with 
density functional theory (DFT) calculations has been extended to describe spin-polarized 
systems. We present calculations on nanocontacts made of Ni as a prototypical magnetic 
material. We find that atomic disorder in the contact region strongly reduces the a pri- 
ori high spin-polarization of the conductance leading to rather moderate values of the 
so-called ballistic magnetoresistance (BMR). On the other hand, we show that the ad- 
sorption of oxygen in the contact region could strongly enhance the spin-polarization of 
the conduction electrons and thus BMR by eliminating the spin-unpolarized s-channel. 
Finally, we show that short atomic Pt chains suspended between the tips of a nanocontact 
are magnetic in contrast to bulk Pt. However, this emergent nanoscale magnetism barely 
affects the overall conductance of the nanocontact making it thus difficult to demonstrate 
by simple conductance measurements. In conclusion, we find that spin-transport through 
atomic-scale conductors is quite sensitive to the actual atomic structure as well as to the 
chemical composition of the conductor. This presents both, opportunities and challenges 
for the realization of future nanoscale spintronics devices. 
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1 Introduction 



The invention of the integrated circuit (IC) in 1959 by J. Kilby (Nobel price in 2000 
together with H. Kroemer and Z. 1. Alferov) triggered the stunning and still ongoing 
development of computer technology which has lead to ever faster, cheaper, and smaller 
computers. An IC is an electronic circuit where all electronic components (transistors, 
capacitors, interconnects) are integrated on a single silicon chip. The successive im- 
provement of the fabrication techniques and the introduction of new materials allowed to 
constantly decrease the sizes of the IC components so that an increasing number of them 
could be integrated on a single IC. This development has lead to increasingly powerful 
and faster computer chips. The other basic ingredient of modern computer systems are 
non-volatile mass storage devices presented in today's computers by hard drives which 
store data in form of magnetic bits on magnetic disks (hard disks) . Here the improvement 
of fabrication techniques, introduction of new materials and new concepts for the the read 
and write mechanism of the hard drives allowed to constantly decrease the minimum area 
needed to record a magnetic bit making them faster and increasing their data capacities 
by several orders of magnitude since the invention in 1979. 

Thus miniaturization has become the leading paradigm of today's computer industry 
and a lot of effort in industrial research and development is dedicated to further decreasing 
the minimum feature size of semiconductor chips or the areas of the magnetic bits on hard 
disks. However, this miniaturization trend cannot go on forever since the atomic scale 
presents an ultimate limit which can never be surpassed and possibly not even reached. 
This is not a kind of hypothetical scenario but rather relevant for the development of 
microprocessor chips in the near future as is impressively demonstrated by today's most 
advanced microprocessors and memory chips whose smallest features (the gate length) 
have already reached a size in the order of 30 nm, and IC transistors with a gate length 
of only lOnm are currently under investigation [Tj. The miniaturization of IC compo- 
nents has in fact reached a level where atomic structure effects like electro -migration 
processes [2] become an issue because they can seriously limit the usability and lifetime 
of ICs. Another serious problem are the so-called subthreshold leakage currents between 
the electrodes of the transistor which increase considerably as the size of the transistors 
shrinks and are responsible for a grand part of the power consumption loss in microchips 
[5]. Similar problems arise when reducing the dimensions of magnetic bits. Moreover, 
when finally reaching the nanoscale, inevitably quantum effects will come into play and 
seriously challenge current silicon based IC technology possibly demanding radically dif- 
ferent approaches. Molecular electronics and spintronics are two promising examples of 
radically new strategies for information processing. 

Molecular electronics aims at employing single molecules as the ultimate electronic 
components for realizing nanoscale electronic circuits. Indeed the use of single molecules 
as rectifiers for electrical current has been proposed as early as 1974 by Aviram and Ratner 
d] . But not until the advent of the scanning probe microscopes was it possible to study 
and manipulate material properties at the molecular or atomic level let alone electrically 
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1 Introduction 

contact individual molecules. The invention of the scanning tunneling microscope (STM) 
by G. Binning and H. Rohrer in 1981 (Nobel price in 1986 together with E. Ruska) |S] 
made it for the first time possible to study (metallic) surfaces and molecules adsorbed on 
them with atomic resolution. By contacting an individual molecule adsorbed on a metal 
surface with the STM tip it is possible to measure the conductance of this molecular 
conductor. This technique was first employed to measure the conductance of individual 
Ceo molecules [SI E] . Since then more conductance measurements of individual molecules 
have been reported in the literature [HI [HI UHl [HI [HI- However, establishing electrical 
contacts with individual molecules and measuring their conductance remains a formidable 
task and care must be taken to assure that a molecule has indeed been contacted [E]. 

Using an STM it is also possible to fabricate atomic-size nanocontacts where two sec- 
tions of a metal wire are connected via a constriction of just a few atoms in diameter which 
thus represent the ultimate electrical conductors with respect to size [T3]- Nanocontacts 
are formed when an STM tip is pressed into the substrate and then slowly retracted 
until an atomic-size neck is formed |14j . Another technique to fabricate very stable and 
reproducible nanocontacts in an efficient and cheap way is the mechanically controllable 
break junction technique (MCBJ)[Tn| where a notched thin wire mounted on a bending 
beam is broken in controlled way. The bending can be fine-controlled by a piczo cle- 
ment allowing a very precise control over the separation of the two sections of the wire. 
Just a third method to fabricate nanocontacts is by electrodcposition [TH [17] where the 
nanocontacts arc clectro-chemically deposited between two macroscopic electrodes. In- 
teresting phenomena for atomic-size nanocontacts have been observed in experiments like 
the formation of monatomic chains suspended between Gold or Platinum contacts, see 

e.g. m- 

Another approach that promises to revolutionize conventional electronics is the field of 
spintronics |19| which aims to combine the traditionally separated fields of magnetic infor- 
mation storage and semiconductor electronics in order to build more powerful electronic 
devices that exploit the electron spin in addition to the electron charge. The example 
that best illustrates the spintronics philosophy is the so-called Magnetic Random Access 
Memory (MRAM) [20] which combines the advantages of conventional magnetic data 
storage (hard drives) and conventional electronic random excess memory (RAM) into a 
single device in order to achieve at the same time non-volatile memory cells that are as 
fast as conventional RAM cells. 

An essential ingredient for spintronics is the generation of spin-polarized electron cur- 
rents which is typically accomplished by passing the electrical current through a ferro- 
magnetic metal. The other basic ingredient are spin- valves which are devices that change 
their resistivity depending on the polarization of the spin-current and thus allow to detect 
the spin-polarization of an electrical current. An example for an actual spin- valve device 
are the giant magneto-resistance (GMR) devices consisting of alternating magnetic and 
non-magnetic metal multilayers that display a strong sensitivity of the electrical current 
to the relative orientation of the magnetizations of the magnetic layers [H] . Soon after 
its discovery in 1988, the GMR effect was exploited to improve the sensitivity of read 
heads in hard drives which before had been based on simple magnetic induction or the 
much smaller anisotropic magneto-resistance effect displayed by bulk metals. This in turn 
allowed to decrease the size of the magnetic bits and thus to increase the data storage 
density of the hard disks dramatically. Consequently, GMR read heads can now be found 
in the hard drives of every modern computer. 



10 



Magnetic tunnel junctions (MTJs) arc similar to GMR devices but feature an insu- 
lating layer instead of the non-magnetic metal layer separating the ferromagnetic metal 
layers which presents a tunnel barrier for the electrons flowing between the ferromagnetic 
layers. MTJs display a MR known as tunneling magneto-resistance (TMR) which was 
first demonstrated by JuUiere [22] ■ The TMR spin- valve is a crucial ingredient for an 
efficient realization of the above described MRAM device [23] . In the earlier experiments 
TMR was found to be much smaller than GMR, but very recently new material combi- 
nations (Fe-MgO-Fe) for the MTJs motivated by theoretical studies [211 12S] have lead to 
a dramatic increase of the TMR which now actually exceeds GMR values [23 123] • 

GMR and TMR spin-valves are nanoscale devices in the sense that the thicknesses of 
the layers making up the spin-valve structures are on the nanometer scale, and can be 
made as thin as 5A with modern growth methods, so that electron transport is governed 
by quantum effects. In fact, these devices actually only work because of quantum ef- 
fects, i.e. the coherent transmission of a spin-polarized current through the non-magnetic 
layer. Moreover, when finally shrinking the other two dimensions of a spin-valve to the 
nanoscale, quantum and atomic structure effects will become even more important. Thus 
it is of fundamental importance to gain a thorough understanding of the interplay between 
magnetism and electrical conduction at the atomic scale. Nanocontacts and nanowires 
made from ferromagnetic metals allow to study this interplay between magnetism and 
electrical conduction or in more fundamental terms the interplay between electron spin 
and charge flow in the smallest possible magnetic conductors. 

An important question is whether GMR or TMR effects survive when the other two 
device dimensions are scaled down to the nanoscale, or whether other magneto-resistance 
(MR) effects emerge at the nanoscale which could be exploited for the realization of 
nanoscale spintronics devices. Therefore, measuring the MR of ferromagnetic nanocon- 
tacts has recently attracted a lot of interest [27 l [28 l [29 l [30 l [n ] [32 ] [33 ] [M l [M l [36 l [37 l [^ 
Indeed, some groups have found a huge MR (exceeding even the GMR effect) for Ni 
nanocontacts which was coined ballistic magneto-resistance (BMR) for its supposed ori- 
gin in the ballistic scattering of spin polarized electrons on a sharp domain wall (DW) 
which should form at the atomic neck of the nano contact [55 ] HP ] 141 j. However, the pos- 
sibility of huge BMR in ferromagnetic nanocontacts has been a controversial topic since 
its discovery, and is one of the principal topics of this thesis (See Ch. |S]and Ch. |S]). 

Integrating the fields of molecular electronics with that of spintronics is another promis- 
ing approach for realizing nanoscale spintronics devices because of the expected very long 
spin-decoherence times in organic molecules as compared to bulk metals and semiconduc- 
tors. Strong indications on the possibility of integrating the two fields come from a few 
recent experiments indicating e.g. a long spin-flip scattering length for carbon nanotubes 
[42], spin injection from strongly spin polarized materials into carbon nanotubes [43], and 
spin transport through organic molecules |44| . Since ferromagnetic nanocontacts arc a 
basic ingredient for molecular spintronics for injecting spin-polarized currents into the 
molecules, it is also in this context important to gain a solid understanding of electrical 
conduction and magnetism in nanocontacts. 

In this thesis electron transport through magnetic nanocontacts and nanowires is in- 
vestigated theoretically. The starting point for the theoretical investigation of electron 
transport through nanocontacts and nanowires is the Landauer formalism, which is in- 
troduced in Ch. 121 The Landauer approach assumes that electron transport through 
nanostructurcs is phase coherent, i.e. decohcrence by phase-breaking scattering processes 
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1 Introduction 

is neglected. This turns out to be a reasonable assumption at low temperatures and for 
small bias voltages. Indeed, the Landauer approach has been successfully applied for 
studying electrical transport through metallic nanocontacts, so that now we have a good 
understanding of atomic scale conductors [ISl US] . 

In order to predict the transport properties of atomic scale nanocontacts and nanowires 
it is important to have a realistic description of the electronic and magnetic structure of 
the nanocontact taking into account its actual atomic structure. This can be achieved 
most conveniently by ab initio electronic structure methods based on localized atomic 
orbitals like e.g. GAUSSIAN [17] or SIESTA gS]. Ch. E] shows how ab initio quantum 
transport calculations based on density functional theory (DFT) [49] are implemented 
in the ALACANT package [50l [51] which as a part of this thesis was extended in order 
to describe spin transport in magnetic systems and to incorporate one-dimensional leads 
calculated from first principles in addition to the semi-empirical Bethe lattice electrodes 
|52j . Other DFT based quantum transport methods are developed by various groups 
around the world [531 [M]. In Ch. [U the Landauer formalism is generalized to describe 
transport of spin-polarized electrons (or in short "spin transport" ) in magnetic nanostruc- 
tures. In order to understand basic aspects of spin transport in nanoscopic conductors 
like the scattering of spin-polarized electrons on domain walls as well as the formation of 
domain walls in nanoscopic conductors, the spin-resolved transport formalism is applied 
to simplified models of magnetic materials. 

In Ch. [5| spin transport through Ni nanocontacts is investigated theoretically with ab 
initio quantum transport calculations using the afore mentioned ALACANT package. In 
order to asses the above discussed possibility of huge BMR in Ni nanocontacts, we cal- 
culate the magneto resistance due to the formation of a DW at the atomic neck of the 
nanocontact for different contact geometries. We find that BMR of pure Ni nanocon- 
tacts is rather moderate [55], i.e. much smaller than the famous GMR effect. This is 
contrary to the claims of huge BMR found in the first experiments on Ni nanocontacts 
but in agreement with some recent experiments measuring very clean samples under very 
controlled conditions |371 [35] . 

In Ch. [H we study the electronic and magnetic structure, and transport properties 
of atomic NiO chains, both ideal infinite ones and short ones suspended between the tip 
atoms of Ni nanocontacts. It is found that the presence of a single oxygen atom between 
the tip atoms of a Ni nanocontact increases the spin-polarization of the conduction elec- 
trons dramatically, converting the nanocontact into an almost perfect half-metallic con- 
ductor, and leading to huge values of the BMR [5B]. It is also discussed to what extent 
these results could explain the huge MR of Ni nanocontacts obtained in the experiments 
mentioned above. 

In Ch. [3 the electronic structure and transport properties of magnetic Pt nanowires is 
studied. It had been shown before that atomic Pt nanowires can actually become magnetic 
in contrast to bulk Pt due to the lower coordination of the Pt atoms in the atomic chain 
compared to the bulk [571 [S5]. We reproduce this result for infinite chains, and find that 
also short Pt chains suspended between the tips nanocontacts become magnetic. However, 
the overall conductance of the Pt nanocontact is barely affected by the magnetism of the 
chain, so that simple conductance measurements of Pt nanocontacts cannot probe the 
magnetism |59| . 

Finally, Ch. [5]concludes this thesis with a general discussion of the obtained results, and 
their significance in the broader context of spin transport and spintronics in nanoscopic 
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systems. Moreover, we will point out open questions and give an outlook on possible 
future lines of work. 
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2 Quantum theory of electron transport 



The description of electrical conduction through a nanoscopic conductor like a molecule 
bridging the tips of two metal electrodes or a nanoscopic constriction in a metal wire as 
shown in Fig. 12. H is a challenging problem. In systems of such small size the dimension of 
the conductor becomes comparable to the Fermi wavelength of the conduction electrons, 
so that the transport properties of the conductor are governed by quantization effects 
demanding a full quantum treatment of the transport process. Moreover, for molccular- 
and atomic-size conductors the actual atomic structure of the conductor has a strong 
effect on the electronic structure and transport properties. 

Our starting point for a quantum description of electrical conduction is the Landauer 
formalism [SUl [HI HH [Ml IM] which is introduced in Sec. 12.11 In Sec. I2.2l wc describe the 
Landauer formalism within the framework of non-equilibrium Green's functions (NEGF) 
formalism. In Sec. I2.3l wc illustrate the derived formalism by applying it to simple model 
systems. Finally, in Sec. 12.41 we discuss the validity of the Landauer approach, pointing 
out its problems and limitations, and indicating ways of improving it as well as alternative 
approaches. 



2.1 Landauer formalism 

In the Landauer formalism electron transport is considered as a scattering process where 
the nanoscopic conductor acts as a quantum mechanical scatterer for the electrons coming 
in from the leads. It is further assumed that the electrons scatter only elastically on the 
nanoscopic sample, i.e. inelastic scattering e.g. by phonons or by other electrons is 
neglected so that transport becomes phase coherent. Thus electron transport through a 
nanoscopic conductor is described in terms of non-interacting quasi particles coming in 
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Figure 2.1: Sketch of typical nanoscopic conductors. Left: A molecule connecting the tips 
of two bulk metal electrodes. Right: A nanoscopic constriction of a metal wire connecting 
two bulk metal electrodes. 



15 



2 Quantum theory of electron transport 




Figure 2.2: Left: Schematic illustration of Landauer formalism (See text for explanation). 
Right: Energy dispersion of propagating modes for electrons moving freely in one direction 
and confined in the other directions. 

Fig. 12.21 shows how the transport problem is modeled in the Landauer formalism: The 
central scattering region (5) containing the nanoscopic conductor is connected via two 
ideal semi-infinite leads of some finite width to two electron reservoirs that are each in 
thermal equilibrium but at different chemical potentials and ^r. The reservoirs are 
assumed to be reflection-less, i.e. incoming electrons are not reflected back to the leads, so 
that the two reservoirs are independent of each other. Due to the finite width of the leads 
the motion of the electrons perpendicular to the direction of the leads is quantized giving 
rise to a finite number of propagating modes or bands {tpn.k} for a given energy E. This 
is illustrated on the right hand side of Fig. 12.21 for the case of a two-dimensional electron 
gas confined in the vertical direction but moving freely in the horizontal direction. In 
that case a propagating wave is given by a plane wave in the z-direction e*'^^ modulated 
by a transverse wavefunction tpnix, y): 

In the more general case of a periodic potential, V{x, y, z + a) — V{x, y, z), generated e.g. 
by the atomic nuclei in the direction of the lead the propagating modes are Bloch waves 

ilj„.k{x,y,z) = Cn-a (fc) 0j;a {x, y, z) e'^"-^ , 

where the sum goes over all unit cells j of the lead and (j)ja{x,y, z) is a localized wave- 
function centered in unit cell j. From the fact that (j)ja{x,y,z + a) = (j)j^ia{x,y, z), it 
follows that the Bloch functions also have the periodicity of the potential, apart from a 
trivial phase factor: 

^pn■,k{x,y,z + a) = '^Ca-n{k) (l)j;a{x,y, z + a) e''^"-' 

= e^'^" Yl ^J-^A^' V' ^) e'""^^-^^ = e^'^'M^, y, z). 
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2.1 Landauer formalism 



We define kn{E) as the wave vector corresponding to the band n at the energy E 
that gives rise to a current in the positive z-direction. The current associated with the 
propagating mode n is given by its group velocity which in turn is given by the derivative 
of the dispersion relation \ for that band: 

,„(fc):=^t;„(fc) = ^^(fc), (2.1) 

where A is the length of the conductor. This is trivial to proof for the case of free electrons 
where the group velocity is directly proportional to the wave vector, Vn{k) = Hk/rrie. For 
Bloch waves the situation is more complicated, and the group velocity can even have the 
opposite sign of the wave vector. A proof can be found e.g. in the book by Ashcroft and 
Mermin [55] . 

Elastic scattering means that an electron with some energy E coming from one of the 
reservoirs will be scattered to some out-moving state with the same energy E of one of 
the leads, so that the process is phase coherent. Thus an incoming wave at some energy 
E on the left lead will give rise to a coherent superposition with outgoing states of the 
same energy E on both leads: 

^n;k„iE)+ E ^""'(^)V^^';-fc„,(B)+ E (^)^n" (B) , (2-2) 

n'eNL n"eNR 

where r„„/ (E) is the probability amplitude for an incoming electron on mode n at en- 
ergy E to be reflected into the outgoing mode n' of the left lead, and tnn"{E) is the 
amplitude for the electron to be transmitted into the mode n" of the right lead. Thus 
an incoming electron on mode n of the left lead will be transmitted with a probability of 
J2n" W^nn" {E)\\'^ to thc right lead giving rise to a current density in that lead of magnitude 

f„{E):^ E \\tnn"{E)rjn"{kn"{E)). (2.3) 

The left electron reservoir injects electrons into thc right-moving modes of the left 
lead up to the chemical potential fj,L- Thus the transmission of electrons through the 
nanoscopic conductor gives rise to a current to the right in thc right electrode. 

/t=E/ dkf„iE,,{k))= E r dEV^,{E)\\t,,^,{EWj,Akn'{E)) 

In the second step we have converted the integral over the wave vectors into an integral 
over the energy by making use of the density of states (DOS) 'D^{E) projected onto band 
n of the right lead. For one-dimensional systems the DOS is given by the inverse derivative 
of thc dispersion relation of the band, I?„(i?) = ^^"i so that it cancels exactly with 
thc group velocity of thc band: 

It= E h dE\\t^r,{EW^-Y. dET^iE), (2.4) 

where we have defined thc transmission per conduction channel Tn{E) as 

T„(S)= E \Kn'{E)f. (2.5) 
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2 Quantum theory of electron transport 



On the other hand the right electron reservoir injects electrons into the left-moving 
modes of the right lead up to the chemical potential /ir and the transmission of electrons 
through the device region gives rise to a left-directed current in the left lead: 



= rdET:SE), (2.6) 



e 

where T,' , {E) now is the transmission probability of channel n' of the right lead 



T'AE)^ J2 K'niEW, (2.7) 

and t'^,^[E) is the transmission amplitude for a mode n' of the right lead to be transmitted 
into mode n of the left lead. Because of time inversion symmetry the amplitude t'^,^[E) 
is the same as the the amplitude tnn'{E) for the transmission from the left to the right 
electrode apart from a trivial phase factor. Hence the total transmission probability from 
the left to the right lead T{E) Tn{E) is equal to the total transmission probability 

from the right to the left lead T'{E) := J2„' K' iE)- 

Y.T.n{E)= J2 \Kn'{E)f= J2 \\Cn{E)f = Tn'iE). (2.8) 

Furthermore the summed reflection probability R'{E) for all electrons injected from the 
right reservoir at some energy E is R'{E) = Nr — T'{E) = Nr — T{E). Therefore the 
current composed of backscattered electrons (originating from the right reservoir) and 
transmitted electrons (originating from the left reservoir) cancels exactly the current of 
the incoming electrons coming in from the right electron reservoir. Analogously the same 
holds true for the left lead. Thus the net current at some energy E is zero when the 
current there is current injection from both electrodes at that energy. 

Thus assuming a positive bias voltage V, so that = + el^ > /xj^, only electrons 
above fXR give a net contribution to the total current. Since for energies above the 
only contribution to the net current through the right lead is the transmission current It 
of electrons coming from the left, the total current for a given bias voltage V is given by 
the famous Landauer formula: 

= J E r dETr^iE). (2.9) 

Taking the derivative with respect to the bias voltage one obtains the corresponding 
conductance: 

G(^) = 5^ - X E ^"(^^)' (2-10) 

2 

where ^ is half the fundamental conductance quantum Gq. The spin-degree of freedom 
of the electrons is contained in the index n for the channels. Assuming a spin-degenerate 
system the transmissions for up- and down channels are equal, so that a factor of two 
appears when summing transmissions over the spin-degree of freedom, and one obtains 
the usual Landauer formula with Gq as the proportionality constant. Here however, we 
are interested in magnetic systems, so the transmissions arc spin-dependent. 
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2.2 Non-equilibrium Green's function formalism 



The transmission amplitudes tnm{E) define an in general non-quadratic matrix t{E) := 
{tnm{E))- The square of this matrix defines a (quadratic) hermitian matrix called the 
transmission matrix: 



The channel transmissions Tn{E) = Tliu'eNn, ^^"^ '^'^'^ J^^t the diagonal ele- 

ments of this transmission matrix, and summing up over all channel transmissions in the 
Landauer formula now corresponds to taking the trace of the transmission matrix: 



The transmission matrix is the central quantity in the Landauer formalism since it 
allows to calculate the electrical conductance and current-voltage characteristics of a 
nanoscopic conductor. Depending on the actual system a variety of methods exists for 
obtaining this quantity. For our purpose of describing transport through atomic- and 
molecular-size conductors the NEGF described in the next section is the most appro- 
priate approach since it can be combined in a straight-forward manner with ah initio 
electronic structure methods like density functional theory (DFT) or the Hartree-Fock 
approximation (HFA) as implemented in standard quantum chemistry codes employing 
atomic orbital basis sets. 

2.2 Non-equilibrium Green's function formalism 

In this section we will restate the Landauer approach to quantum transport in the lan- 
guage of one-body Green's functions [Ml [M|- To this end we will first introduce the 
Hamiltonian and overlap matrices for the transport problem in a basis set of localized 
atomic orbitals. We will then derive a Green's function for the finite scattering region con- 
nected on both sides to semi-infinite leads. From the Green's function (GF) one can then 
obtain the (reduced) density matrix and the transmission matrix. Here we will mainly 
follow the arguments presented by Paulsson in his introductory paper on the NEGF [S7] 
but generalize them to non-orthogonal basis sets (NOBS) as commonly employed in quan- 
tum chemistry packages. A few other derivations of the Landauer formalism within the 
NEGF framework taking into account non-orthogonality of basis sets can be found in the 
recent literature [551 

2.2.1 Hamiltonian and overlap 

We divide the system into 3 parts as shown in Fig. 12.31 The left lead (L), the right lead 
(R), and the intermediate region called device (D) containing the central scattering region 
(S). This scattering region can be given by e.g. a molecule coupled to metallic contacts, 
or simply a nanoscopic constriction as indicated in the figure. 

We assume that the leads arc only coupled to the scattering region but not to each 
other. Thus the device region must be chosen sufficiently large for that to be true. The 




(2.11) 




(2.12) 
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r 




Device D 



Figure 2.3: Sketch of the transport prob- 
lem. L: Left (bulk) lead. D: Device. R: 
right (bulk) Lead. S: Central scatter- 
ing region containing nanoconstriction 
or contacted molecule. See text for fur- 
ther explanations. 



Hamiltonian H describing the system is then given by the matrix 



H 




(2.13) 



Furthermore we assume a non-orthogonal localized basis set. Assuming again no overlap 
between atomic orbitals in different leads the overlap of the atomic-orbitals of the system 
is given by the following overlap matrix: 



(2.14) 



As indicated in Fig. I2.3l we subdivide the leads into unit cells (UCs) which must be chosen 
sufficiently large so that the coupling between non-neighboring unit cells can be neglected. 
Thus in general a UC consists of several primitive unit cells (PUCs). The Hamiltonian 
matrix H^, of the left lead can be subdivided into sub-matrices in the following manner: 



/ 



H 



L — 



-2,-2 ■n._2,-l 

\ ■ ■ ■ H_i^_2 H_i^_i 



Hi Ho Hi 

Hi 



^ ' ' ^ - : / HI Ho 

\0 Hi Ho/ 

Analogously the Hamiltonian of the right lead is given by the following matrix: 



(2.15) 



H, 



^ Hi^i Hi^2 

H2,l H2,2 



/Ho Hi 
Hi Ho Hi 

t M„ Hi 



Hi Ho 



Vo 



(2.16) 
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In a similar way, the overlap inside the leads is given by the matrices 



S_2,-2 S_2,-l 
\ • • • S_i^_2 S_i _i 



Vo 



sl So Sl 

sl So Sl 
sl So/ 



(2.17) 
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and 



/ Si,i 

§2,1 



Sl^2 
§2,2 



Si 
So 
Si 



Si 
So 



Si 



Vo 



7 



(2.18) 



Furthermore the unit cell of each lead that is immediately connected to the scattering 
region (unit cell "I" for the left and unit cell "r" for the right lead) is included into the 
device part of the system.: 







H/.s 




= 




Hs 






^ Or,; 




H, 



and Sd 




(2.19) 



In the next subsection wc will show how to calculate the GF for the device part of the 
system as defined by the above Hamilton and overlap matrices. 



2.2.2 Green's function for the open system 

The one-body GF operator G{E) of a system is defined as the solution to the generalized 
inhomogeneous Schrodinger equation (see e.g the book by E. N. Economou, Ref. [55]): 

(z - i?)G(z) = i, (2.20) 

where H is an (effective) one-body Hamiltonian, and z is a complex number. When z 
does not coincide with the eigenvalues of the Hamiltonian H, the GF operator has the 
following formal solution: 

G{z) = {z ~ H)-^ for z ^ ek. (2.21) 

Obviously, for z = ek the GF operator has a pole and is thus not well defined. In this 
case one can define two GFs which are both solutions to eq. (|2.20p by a limiting process. 
1) The retarded GF is defined as: 

G{E) \im{E + irj- H)-'^, (2.22) 
jj— »o 

and 2) the advanced GF is defined as the hermitian conjugate of the retarded GF: 

G^E) := lim{E - irj - Hy\ (2.23) 

E is a, real number which can be any of the energy eigenvalues of the Hamiltonian 
H. When E does not coincide with an eigenvalue of H the two GFs reduce to the GF 
operator defined in eq. (|2.2ip . In the basis of eigenstates of the Hamiltonian the GF 
operator becomes diagonal: 

G(z) = V with H\k) = ek \k). (2.24) 

fe " 

The GF yields the complete information of a one-body system. For example, eq. (|2.24p 
makes clear that the poles of the GF along the real axis represent the eigenvalues of the 
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Hamiltonian H. Thus by plotting G{E) along the real axis one can find all the eigenvalues 
in a certain energy range. 

A very important quantity is the density of states (DOS) T>{E). The DOS can be 
calculated from the trace of the imaginary part of the GF along the real axis: 

Im Ti\G(E)] = y lim Im ^- = 

k 

k k 

Another important quantity of the GF formalism is the spectral density A{z) which is 
defined as the difference between the retarded and the advanced GF: 

A{E):^i{G{E)-G\E)). (2.26) 

It is straight forward to show that the eigenstate representation of the spectral function 
on the real axis is given by: 

A{E) ■.^Y.^iE - ek)\k){k\, (2.27) 

k 

so that the trace of the spectral function directly gives the DOS. The spectral function 
can thus be seen as a generalized DOS. 

The grand advantage of the GF formalism is that it allows one to calculate all properties 
of a one-body system without having to calculate the eigenstates of the Hamiltonian 
explicitly. Instead the GF can be calculated in any basis set by a matrix inversion for any 
value z, eq. (j2.2ip . It turns out that in many situations this is more convenient than to 
solve the whole eigenvalue problem. This is the case for the transport problem defined by 
the Hamiltonian and overlap matrices given in the previous subsection, eqs. (|2.13p and 

(ma. 

The matrix GapiE) := (^a\G{E)\(3) corresponding to the GF operator in a non- 
orthogonal basis set { | a) } is given by: 

{ES- n)S-^G{E) = S, (2.28) 

where S is the overlap matrix of the NOBS: Sap = (a | /3). This, however, is a somewhat 
inconvenient definition for the calculation of the_^GF since it involves the inversion of the 
S matrix. Instead, we define a new GF matrix G{E) by 

G{E) := S-^G{E)S-^. (2.29) 

This GF can be calculated more conveniently from the much simpler equation 

(ES-n)G{E) = 1. (2.30) 

For convenience we introduce the following energy-dependent matrix: 

Hi— _ESi Hlo^ESlo 
W(£;) :=H-^;S= I Hdl-ESdl Hd^ESd Udr- ESdr |. (2.31) 

HijD — E Srd H/j — E Hr 
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1 






^ 








Thus we obtain for the GF for the transport problem the following matrix equation 
which defines a system of equations for each sub-matrix of the total GF. 

Wl{E) WLDiE) WLRiE) 
Wdl{E) Wd{E) WdbXE) 
Wrl{E) Wrd{E) Wr{E) 

( Gl{E) Gld{E) Glr{E) 
X Gdl{E) Gd{E) Gj^rIE) I = I 1 I , (2.32) 
\ GrlIE) GRDiE) GRiE) 

This matrix equation can be solved for each of the matrix elements (see App. [B]). For 
the device part of the GF Gjj we obtain: 

Gd{E) = {ESd Hd ^l{E) - ^R{E)r\ (2.33) 

where we have introduced the self-energies of the leads Si and which describe the 
influence of the leads on the electronic structure of the device. They can be calculated 
from the GFs of the isolated left and right semi-infinite lead, ^l{E) = (ESl — H^)^-'^ 
and gR{E) ^ {ESr — H^)^^, respectively: 

Si(S) = WDL{E)gL{E)WLD{E) (2.34) 
Sfl(S) = WDR{E)gRiE)WRD{E) (2.35) 

For later use, we also define the so-called coupling matrices: 
fL{E) := z (Si(i?)-S[(i?) 

= iWDL{E){gL{E)~si{E))WLD{E), (2.36) 
TRiE) := t (^LiE)-±i{E)] 



= iWDRiE)igRiE)^g^j,{E))WRDiE). (2.37) 

Thus we have expressed the GF of the device region in terms of the device Hamiltonian 
and the Green's functions of the isolated leads. The term H/j + Sl(-E') + 'Er{E) can be 
interpreted as an effective Hamiltonian for the device region, and its energy dependence 
stems from the fact, that the lifetime of an electron in the device region is now finite due 
to the coupling to the leads. The technique of calculating the GF in parts by dividing 
the underlying Hilbcrt space into subspaccs is called partitioning technique and is shown 
in more detail in App. [Bl where we also list the other matrix elements of the GF. 

Since the coupling Hdl between the left lead and the device is only due to coupling 
between the last unit cell (—1) of left lead and the unit cell / included in the device 
region, only the surface GF of the left lead, i.e. the GF projected into the unit cell — 1, 
is necessary for calculating the self-energy S/,. The same holds true for the right lead in 
an analogous manner. Furthermore, S^, is different from zero only in the Z-region of the 
device, and analogously Hr in the r-region: 




S := Si + Sfl = Os,; Os Os,,. , (2.38) 
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where the non-zero matrix elements and of S can be expressed in terms of the 
surface GFs of the left and right lead, g-i,-i and gi,i: 

±i{E) = (H\-ES\)s-i^-iiE)(Ui~ESi), (2.39) 
flriE) = (H,-ESi)g+i^+i{E){H\-ES\). (2.40) 

The self-energies can be calculated iteratively by Dyson equations, as shown in App. 

El 

%{E) - {H\-ES\){ES„-Ho-i:i{E))-^{Ui-ESi), (2.41) 

f^riE) - {lli-ESi){ES„-Uo-^r{E)r^{H\-ES\). (2.42) 

Now we have achieved a description of the electronic structure of the device region in 
terms of one-body GFs which take into account the effect of the coupling of the leads 
to the device region. In the next section we will show how to obtain the (reduced) 
density matrix and subsequently the electron density of the device region from this GF. 
Thereafter we demonstrate, how to calculate the transmission matrix from the GF of the 
device region and the self-energies of the leads. 

2.2.3 Calculation of density matrix and electron number at equilibrium 

The reduced density matrix of first order (in quantum chemistry often called charge 
density matrix) is obtained by tracing out all but one of the one-particle subspaces from 
the many-body density matrix p: 

P = Nx ^ (n,n2, . . . ,njv|/5|n', n2, . . . , nAr)|ri)(ri'|, (2.43) 

n2,...,njv 

where {|?^)} is a set of one-body states. 

At zero temperature the system will be in its (many-body) ground state so that the 
density matrix becomes 

pv[, = |*)(*|. (2.44) 

It is straight forward to show that in this case the reduced density matrix can be expressed 
as 

P = J2\k){^\clci\^){l\. (2.45) 

k.l 

One can obtain a lot of information from the reduced density matrix about a many-body 
system in spite of the fact that it is actually a one-body operator. For example, one can 
calculate the expectation value of any one-body observable from the trace of the product 
of the reduced density matrix and the observable: 

(i)=Tr[Pi]. (2.46) 

The electron density is given by the diagonal elements of the reduced density matrix in 
the real space representation: 

n{r) = {r\P\r)- (2.47) 
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Finally, the number of electrons is the trace of the density matrix. 

In the Landauer approach we consider the electrons as a system of non-interacting 
quasi-particles, i.e. the Coulomb interaction between the electrons is only taken into 
account on the mean-field level. In the case of a system of non-interacting particles the 
many-body ground state of the system j^*) is given by a single Slater determinant: 

|*)= n 4|0), (2.48) 

fe,efc</i 

where are the eigenstates of the one-body Hamiltonian, |o) is the vacuum ground 
state, and fi the chemical potential, so that the Slater determinant consists of all states 
with energies less than or equal to the chemical potential. In this case the reduced density 
matrix is diagonal in the one-body eigenstates \ky. 

P = J2nk\k){k\, (2.49) 
fc 

where the occupation number of the eigenstate is given by the Fermi distribution 
function, Uk = /(efe — m)i cind. thus is either one for states below the chemical potential 
or zero for states above the chemical potential at zero temperature. Since we will make 
use of the reduced density matrix only but not of the full many-body density matrix we 
will refer to the reduced density matrix for the sake of simplicity from here on simply as 
the density matrix, unless otherwise stated. 

Using the eigenstate representation of the GF, cq. (|2.24p . it is straight forward to show 
that the density matrix can actually be obtained from the GF of the system by integrating 
the imaginary part of the GF up to the chemical potential fj, of the system: 

1 /-M } IfcVfcl 
Im / dEG{E) = lim / dElmY^ L^^J 

dE \im-— ^- ?|fc)(fc|-y / dES{E~ek)\k){k\ 

= J2fi^k-f^)\k){k\=P. (2.50) 

k 

Subsequently, the (standard) density matrix in a NOBS is obtained by integra- 

tion of the standard non-or thogonal GF matrix, Gap{E) := {a\G{E)\l3): 

Pap ■■= {a\P\f3) = --Im r dEG^p{E). (2.51) 

^ J — oo 

Analogously to the GF matrix G we can define a new density matrix by 

P = S^PS\ (2.52) 
which is thus obtained by integration of the non-standard GF G(£'): 

P = Im / dEG{E). (2.53) 
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One can obtain the number of electrons from the trace of the density matrix. However, 
the trace of an operator has to be taken in an orthogonal basis set, but P has been defined 
in terms of the density matrix in a NOBS The standard procedure for orthogonal- 

izing atomic orbitals is the so-called symmetric orthogonalization |70j . A density matrix 
P in a NOBS will be transformed to the density matrix P-'- in the orthogonalized basis 
according to: 

p± ^ s-l/2 p s-1/2 = ^ + 1/2 p s + (2.54) 

In the last step we have used the definition of the non-standard density matrix P to 
obtain the transformation to the orthogonalized density matrix. 

Applying the above transformation for the density matrix P we find for the number of 
electrons the following expression: 

iVe = Tr[Si/2 p ^1/2] ^ Tr[P S] = Trpc S^,] + Tr[VDL Sld] + Tr[Pz,fl Srd] 

+ Tr[PL Si] + Tr[VLD ^dl] + Tr[Pfl S^] + Tr[Pflz? ^dr] (2.55) 

where we have first exploited the invariance of the trace under commutation of matrices 
and then applied the above explained division into sub-matrices of the density matrix and 
the expression (|2.14p for the overlap matrix. In analogy to the MuUiken analysis we can 
identify the number of electrons in the device, Nd^ in the left lead, N^, and in the right 
lead Nr as follows: 

Nd = Tr[Pz3Sz3] + Tr[PDLSLi,] + Tr[Pz5fl^Sfl.^Z5] 
Nl = TrpL Sl] + Tr[Pi,5 Sdl] 

Nr = Tr[Pfi SflJ + TrpKB Sda] (2.56) 

Thus due to the overlap of the device orbitals with the lead orbitals one also has to 
calculate the off-diagonal elements DL and DR of the density matrix P in order to 
calculate the number of electrons in the device part. However, computing the off-diagonal 
elements of the density matrix would require calculation of the corresponding off-diagonal 
elements of the GF matrix G{E) in terms of the GF matrix of the device and the self- 
energies (see App. [B]) which are quite tedious expressions. Instead we will require charge 
neutrality only for the scattering region and not for the entire device region. That is the 
reason why we have included one unit cell of each lead into the device region in addition 
to the scattering region. Since the scattering region docs not have an overlap with the 
two leads it is sufficient to calculate the density matrix of the device region to obtain the 
number of electrons inside the scattering region: 

Ns — Tr5[P£i S^)] = 

Po.f3 S(3o. (2.57) 

aes.peD 

Finally, we also have to define an appropriate DOS V and projected DOS (PDOS) Vi 
for the case of non-orthogonal orbitals. In order to be coherent with the above definition 
of the electron numbers in the respective subspaccs of the leads and the device, it is best 
to adapt a similar MuUiken-likc definition for the PDOS: 

V,,{E) = --lm{G{E) S)„„ (2.58) 
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where a denotes some of the atomic orbitals of the system. If we are interested in the 
DOS projected onto some orbital a of the scattering region, it suffices to calculate the 
GF and overlap matrix for the device region instead for the entire system, VaiE) = 
— -lm{GD{E) S£))aQ. And the DOS projected onto the entire scattering region is just 

2.2.4 Non-equilibrium density matrix 

In the previous section wc have assumed that the system is in thermal equilibrium, i.e. 
the chemical potentials of the two electron reservoirs are equal: /i^ = M_r = Although 
in many cases one can estimate the conductance at small bias from the zero-bias conduc- 
tance, i.e. in equilibrium, it would be very interesting to study the interplay between an 
electron current and the electronic structure of the nanoscopic conductor. Thus we have 
to derive an expression for the density matrix when the system is out of equilibrium, i.e. 

First, we calculate the response |V'^) of the entire system to an incoming wave on the 
left lead I?/;*) which we assume is a solution of the isolated left lead, i.e. Hl\iP^) = E\ip'^'). 

H{\i^+) + = + IV-')) (2.59) 

Because of the non-orthogonality of the lead orbitals with the device orbitals we now 
switch to a matrix presentation of the problem as done in the previous section. The 
incoming wave as a solution of the left isolated lead is given as an expansion over the lead 
orbitals: 

/ V4 \ 

= 5] Ch) ^ = Od (2.60) 

"GL \ Or J 

The response wave function on the other hand expands over the whole system: 

^+ = I I (2.61) 

Thus in matrix representation the above Schrodinger equation for the response to an 
incoming wave reads: 

H(iA+-hV'') = E^{^f+ + ^f^) 

I {ESL-iiL)A 

=^ (H - ES)i:+ = {ES - = {ESdl ~ iiDLWi 

\ Or 

-G{E) I [ESdl - ^dlWl ■ (2-62) 
Or J 

Since is a solution for the isolated left lead, 

VLlA - ESl^^I ^ (Hi - ESl)^1 = (2.63) 
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we obtain finally: 





/ 


Gld{E) \ 


y 




Gd{E) 




\ 


Grd{E) j 



rf+ = -G{E) (ESdl - JioM 

J V Or 

The left electron reservoir fills the incoming electron waves \ipL,n{k)) on the left lead 
up to the chemical potential /i^ of that reservoir which result in wavefunctions expanded 
over the whole system |\E'„(fc)) = \ip^{k)') + \ipL,nik)). Thus the density matrix due to 
injection of electrons from the left reservoir is given by: 

^^''^ = E /(ei."(fc)-ML)|*„(fc))(*n(fc)| 



J2 I dE6{E - eL.nik))f{E - AiL)|*„(fc)))(*„(A:)| 



I 



dEfiE~fiL) J2 |*n(A;)))5(S-eL,„(A;))(^n(fc)| (2-65) 

n£NL,k 



In matrix representation this density matrix is given as 

= = E /(ei.-W--"i)^"a'^„a'(fc)*„/3'(fc)*V/3 (2-66) 

n£NL,k 

Thus the non-standard matrix representation is 

p{L) = S-1PS-1 = ^n{k)S{E-eL,nm[^n{k)]^. (2.67) 

n<£NL,k 

From eq. (12.6411 we see that the device part is given by 

P^''^ = fdEf{E-fiL) ^DAk)S{E~eL,n{k))[^DAk)]^ 

= j dEf{E-IJL)GDiE)WDLiE) 

J2 i^LAk) S{E - eLAk)) ii^LAk)]^ Wld{E) &AE) 

n&NL,k 

" V ' 

5i(£;)/27r=i(gi,(i5)-i+(£;))/27r 

= ^ I dEfjE - f,L)GAE) {i:L{E) - f:[iE)) GUE), (2.68) 

where in the 2nd step we have made use of the spectral density of the isolated left lead, 

aAE)^^igL{E)-gim■ 

For the injection of electrons from the right we obtain an analogous expression. Sum- 
ming up both contributions we obtain the density matrix out of equilibrium: 

= p(^^) + p(^«) = ± 1 dEGD{E)[f{E ^iL)TAE) + f{E ^^n)rn{E)]GUE) 

(2.69) 
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Finally, we define the non-equilibrium GF matrix G^(£') of the device as: 

G<{E) := tGD{E)[f{E - fiL)rL{E) + f{E - fiR)rniE)]GUE) (2.70) 
Then the non-equilibrium density matrix can be written as: 



p„eg ^ _± I ^EG<{E) (2.71) 

2tt J 

In equilibrium, i.e. when /i^ = = fi, the GF matrix G^{E) reduces to the imaginary 
part of the retarded GF matrix: 

G<{E) = f{E~ ^i)iGD{E)[fL{E)+fR{E)]&JJ{E) 

= .f{E-^,)GD{E)[{GD{E))-' - {&^{E))-^]gUe) 

= f{E-ti)[G^^{E)~Gn{E)] = -2if{E^li)lin[GD{E)l (2.72) 

and one arrives again at the expression for the equilibrium density matrix derived before: 

P = -^y" dEf{E-^L)lm[GD{E)]. (2.73) 

The same reasoning shows that for an out of equilibrium situation the density matrix 
can be calculated using the much simpler expression (|2.73p for the energy integration up 
to the lowest of the two chemical potentials, e.g. for > we have: 

Pd' = r dElm[GD{E)] - ^ ^ dEG<{E). (2.74) 

J-oo 27r J^^ 

2.2.5 Current and transmission 

Now we derive an expression for the electrical current through the device region in terms 
of the GFs by calculating the time derivative of the electron charge in the device region. 
The temporal change of the number of electrons inside the device region is equal to the 
sum of the currents from the left and from the right lead: 



dF 



= Il+ Ir, (2.75) 



where we have taken a current as positive when it increases the number of electrons inside 
the device region, i.e. when the current flows from the leads to the device. Typically we 
are interested in steady-state situations, so that the electron charge of the device will be 
conserved, meaning that II = —Ir- Then the net current / through the device is given 
by / = /l = -Ir. 

We will now derive an expression for the total current into the device ^^f- = Il+Ir by 
taking the time derivative of the density matrix P, and divide the obtained expression for 
the total current II + Ir into the device into its individual contributions II and Ir. We 
start from the expression (|2.67p in order to calculate the net current through the device 
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due to the injection of electrons from the left electrode: 



(9N 



dt 



dt 



(2.76) 



xTr 



dt 



For the first term of the trace we find: 

1 



Tr 



D 



dt 



-Tr 



D 



H*„(fc)[*„(fc)]t 



1 



iiDL[^n{k)]L[^n{k)]l + Hfl (fc)] ^ [*„ (fc)] + Hz^fl (fc)] [*„ (fc)] 



(2.77) 



where we have made use of the invariance of the trace under cyclic permutations and the 

fact that ^„(fc) is an eigenstate of the system, so that ih S '^'^q^^'^ = H \i?„(fe). 
Analogously, we find for the second term of the trace (|2.76p : 



Tr, 



dt 



in 



= -- rTrz, 
in 



*„(fc)[^„(fc)]^H 



[*„(fc)]z5[*„(fc)]LHLI, + [*„(fc)]z5[*„(fc)]T,Hz5 + [^nmoV^ n{k)\\liRD 

(2.78) 



When adding up the two contributions to the trace in (|2.76p . the term involving H^i 
in ()2.77p cancels with the term involving H/^ in (|2.78p . For later convenience we add the 
following term which is equal to zero to the trace: 







It,. 

in 
in 
in 



i;„(fc)s*„(fc)[*„(fc)]t 



E^{k)SDL [^n{k)]L[^n{k)YD 



[^„(fc)]i5[*„(fc)]ii;„(fc)s 



LD 



in 
in 



En{k)SDR [*„(fc)]fl[*„(fc)]T, 
[^n{k)]D[^n{k%Er,{k) 



where the two terms involving Sd have canceled out each other. Summing up all con- 
tributions to the trace in (|2.76p and grouping together on the one hand terms involving 
hopping between the left electrode and the device and on the other hand terms involving 
hopping between the right electrode and the device, we get in total: 



Tr, 



s^i*.(.)|. + $4*,2*5W!!s 



dt 



dt 



= jL.n{k) + jRji{k) 



(2.79) 



where 



jL,n{k) ■■= -rTr ['i'n{k)]lWDL{En{k)) [^n{k)]L - (fc)][ W^z, (£„ (fc) ) [*„ (fc)] , 
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and 



The term jL,n{k) gives the current coming from the left lead, while j.R,„(fc) gives the 
current coming from the right electrode. We proceed by expressing the states |^'„(A:)) in 
terms of the incoming waves |5'Jj(fc)) by the help of expression ()2.64p for the response to 
an incoming wave: 



[^n{k)]^ = GDiE)WDLiE)[Kik)] 
[^n{k)]l = [Kik)]^WLD{E)GUE) 
[*n(fc)]zj = GRD{E)WnL{E)[Km 



= SR (E) Wrd (E) Gd {E)Wdl (E) (fc)] 
[*n(fc)]fl = [^i{k)]^WLDiE)GUE)WDRiE)gUE) 



jR,n{k) = -zTr 
in 



[%m^WLD{E) gUe)x 
Wdr{E) [gR{E)-g^j,{E)} WuDiE) Gom^DLiEWnik) 



(2.80) 
(2.81) 

(2.82) 
(2.83) 



where the energy argument of the GFs E is the band energy of the right lead: E = E^{k). 
Inserting these expressions into the expression for jR^n{k) defined above we get: 

1 



-iTniE) 

&j,{E)fR{E)GD{E)WDL{E)^\{kWn{k)]^WLD{E)] . (2.84) 



Thus summing over all states n resulting from electron injection from the left electrode 
and integrating over energy we get the contribution /j^^' to the current Ir = / 
resulting from electron injection from the left electrode. 



(i) , j(R) 
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/ dE f{E - ^1L)TI\G^^{E)TR{E)GD{E)TL{E) 



(2.85) 



In a similar way we get for the current if^'^ through the left lead resulting from electron 
injection from the right reservoir: 



AR) 



-- I dE f{E - hr,)Tt Gl{E)TR{E)GD{E)TLiE) 



(2.86) 
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In steady-state the contributions to the current due to electron injection from the left 
reservoir on the one hand and from the right reservoir on the other hand are individually 
conserved, i.e. / 



= /(^) and li"^ = -/L'"^' = /(^). Thus in total the net 



current through the device is / = /'•^■' 
famous Landaucr formula for the current expressed with GFs.: 



j(B.) ^ _jiL} 



and we finally obtain the 



dE iUiE) - fR{E)) Tr gUe) Tn{E) Gd{E) Tl{E) 



(2.87) 



The trace in ()2.87|) corresponds to the total transmission function T{E) which is the 
sum over all channel transmissions Tn{E) defined earlier: 



T{E) 



T^{E) = Tv G^j,{E)Tr{E)Gd{E)Tl{E) 



(2.88) 



This is the so-called Caroli expression for the transmission function T{E) named after C. 
Caroli who first derived it for the simple case of a one-dimensional tight-binding chain 



Obviously, the (non-hermitian) expression inside the trace of (|2.88p must have to do 
with the (hermitian) transmission matrix T{E) defined earlier in eq. (|2.1ip . Exploiting 
the invariance of the trace under cyclic permutations we find an equivalent expression for 
the transmission function now involving a hermitian expression that we can identify with 
the transmission matrix: 



T{E) ^ t'j^\e) Gl{E) Tn{E) Gn{E) r'^'{E) 



1/2, 



Now we can also identify the transmission amplitude t{E) = {tmn{E)) with: 



t{E) 



T',^{E)GD{E)T'r{E). 



^1/2/ 



(2.89) 



(2.90) 



The generalized power of a hermitian matrix A is defined as = Udiag(af)U^ where 
the Oi are the eigenvalues of the matrix A and U is the unitary transformation that 
diagonalizes A. 

We note that the dimension of the transmission matrix (|2.1ip obtained in Sec. 12.11 is 
equal to either the number of modes of either the left or the right lead, which in turn is 
equal to the dimension of the unit cell of either lead. But the transmission matrix derived 
here has the dimension of the device subspace which is bigger than the dimensions of the 
unit cells of cither lead. So the transmission matrix (|2.89p derived from the NEGF is 
strictly speaking not equal to the original one (|2.1ip of the Landauer formalism. How- 
ever, both formulations are of course equivalent. In fact, all non-zero matrix elements of 
the transmission matrix in the NEGF correspond to the original Landaucr transmission 
matrix. All other matrix elements of the NEGF transmission matrix arc zero. 



2.3 Application to simple models 

In order to illustrate some aspects of the NEGF formalism derived above, it is applied 
here to some simple model systems. We start with the well-known example of the simple 
tight-binding chain with one orbital per atom. While the standard tight-binding model 
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neglects the overlap between neighboring atomic orbitals, the overlap is included here as 
an additional parameter s in order to illustrate its effect on the electronic structure. Thus 
we have for the matrix elements of the Hamiltonian: 

{eo for m = n, 
t for \m-n\ = 1, (2.91) 
otherwise 

and for the overlap matrix elements: 

{1 for m = n, 
s for |m — n| = 1, (2.92) 
otherwise 

Note, that the representation of the Hamiltonian in operator form in the NOBS is not 
straight forward, since it involves the inversion of the overlap matrix S ~ (Smn) as pointed 
out in App. [Al eq. (|A.2p : 

= ^|m)(S-iHS-i),„„(n|. (2.93) 

mn 

Thus one introduces hopping terms beyond the first neighbor into the Hamiltonian oper- 
ator although the corresponding matrix elements arc zero. Orthogonalizing the basis set 
for example by the Lowdin scheme one would also obtain hopping terms beyond the the 
nearest neighbor hopping. But we will not pursue this any further. Instead we will apply 
the above developed NEGF formalism for a NOBS to the model and study the effect of 
the overlap on the DOS. 

The PDOS for a NOBS is given by eq. (^351) . The PDOS of the tight-binding chain 
is plotted in Fig. 12.41 for various parameters of the nearest neighbor overlap s, while the 
on-site energy and the hopping are kept fixed: eo := and t := —1. We observe, that 
increasing the overlap leads to a displacement towards higher energies and a simultaneous 
broadening of the energy band until at s = 1/2 the band width becomes infinite. This is 
due to a singularity of the energy dispersion at the Brillouin zone boundary ka — ±7r: 

TPfi\ eo + 2icos(fca) r , , -, 1 /o 

E(k) = 7—^ — > oo for ka ±tt and s = -. (2.94) 

^ ^ l + 2scos{ka) 2 ^ ' 

Increasing the overlap s beyond this value, the DOS becomes even negative for some 
values of the energy. This pathological behavior can be traced back to the appearance of 
singularities inside the first Brillouin zone, and relates to the fact that the tight-binding 
model is actually pathological for overlaps s > 1/2. In fact, it is physically impossible 
to have a nearest neighbor overlap of s > 1/2 but no overlap between second nearest 
neighbors. The pathological behavior can be repaired quite easily by including overlap 
beyond the nearest neighbor approximation. This can be seen directly from the energy 
dispersion of the generalized tight-binding model including hopping and overlaps beyond 
the first nearest neighbor: 

. eo + E.2t.cosikna) 
1 + 2s„ cos(fcna) 

where s„ and t„ denote overlap and hopping, respectively, between n-th nearest neighbors. 
Thus in order to avoid unphysical results it is important to include overlap between atomic 
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Figure 2.4: PDOS of tight-binding chain with overlap for different values of the overlap 
parameter s but for fixed on-site energy parameter eo and hopping parameter t. 
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Figure 2.5: Tunneling between two semi-infinite chains. The graphs show the PDOS 
projected onto a tip atom (left) and the transmission (right) for different values of the 
hopping t' between the chains in units of the hopping parameter t of the tiglit-binding 
chains. Overlap between atomic orbitals is completely neglected. 



orbitals up to a sufficiently high degree of neighborhood. This will become very important 
in the context of combining ab initio electronic structure calculations with the NEGF 
formalism as the the basis sets employed in common quantum chemistry packages like 
GAUSSIAN [37] are typically highly non-orthogonal. 

Next, we will have a look at the electronic structure and the transport properties of two 
tunnel-coupled semi-infinite tight-binding chains. In Fig. 12.51 the PDOS projected onto 
one of the tip atoms (bottom left) and the transmission function (bottom right) is shown 
for the model depicted at the top of the figure for various values of the tunnel coupling 
t' between the two chains. We observe, how the PDOS is slowly transformed from the 
DOS of the ID chain with the characteristic van Hoft singularities at the band edges to 
the semicircular "surface" DOS of the semi-infinite chain when the tunneling coupling 
t' is decreased. Simultaneously, the perfect transmission of the ID chain within the 
energy band is reduced considerably when decreasing the tunnel coupling as is expected. 
For t' = O.lt we enter the tunneling regime where the transmission is described by the 
Bardeen formula for timneling T « t' /t'D{eF) [7^ . 

Finally, we have a look at a simplified model of a typical situation in mesoscopic physics 
and molecular electronics: An island (which could be a quantum dot in a mesoscopic 
system or a molecule in molecular electronics) coupled to two semi-infinite leads. The 
island in the simple model is a single site coupled to two semi-infinite tight-binding chains 
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Figure 2.6: Tunneling through a dot connected to two semi-infinite leads. The graphs 
show the PDOS projected onto the central dot (left) and the transmission (right) for 
different values of the coupling t' of the dot to the two leads in units of the hopping 
parameter t of the tight-binding chains. 
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(see the illustration at the top of Fig. 12. 6p . We start again with strong coupling to the 
leads, t' = t, and decrease t'. We observe that the PDOS projected onto the island is 
tranformed from the DOS of the perfect ID chain for {t' = i) to a Lorentz peak for weak 
coupling and also the transmission is now given by a Lorentz peak in that limit. We 
further note that the transmission remains perfect exactly at the on-site energy of the 
island although the coupling has decreased by one order of magnitude. This phenomenon 
is called resonant transmission through a molecular (or atomic) level. 

2.4 Discussion of the Landauer approach 

The Landauer formalism introduced above, assumes that electron transport in nanoscopic 
conductors is phase coherent. Inelastic scattering of electrons e.g. in electron-electron 
or elcctron-phonon scattering processes is neglected. This turns out to be a rather good 
approximation for low temperature and small bias voltages. The reason for this is that for 
small bias voltages the transport properties are principally determined by the conduction 
electrons at the Fermi level. At low temperatures the conduction electrons at the Fermi 
level cannot loose energy because the Fermi sea is completely occupied. On the other 
hand they cannot gain energy because there are no phonons available at low temperatures. 
Therefore electron transport is elastic, and thus phase coherent at low temperature and 
for small bias. 

A strict generalization of the Landauer formula to interacting electron systems (i.e. 
including electron-electron interaction and/or electron-phonon interaction) is given by 
the Meir-Wingreen formula |73j which assumes an interacting device coupled to non- 
interacting leads. The Meir-Wingreen formula actually reduces to the Landauer formula in 
the limit of zero temperature and zero bias, thus demonstrating the above considerations 
strictly. Thus in principle the Landauer formula can be considered a strict result (for low 
temperature and small bias) under the condition that the GF of the device is the true 
interacting GF (see e.g. the book by Mahan [74] for an introduction to many-body GFs). 
However, the GF employed here is the one-body GF describing non-interacting electrons. 
Thus the electron-electron interaction can only be taken into account on a mean-field 
level. The task is now to find a one-body GF that approximates the true interacting 
many-body GF sufficiently well. In Ch. Owe show how to combine the NEGF formalism 
with ab-initio electronic structure calculations on the level of DFT, i.e. the one-body GF 
is calculated from the Kohn-Sham Hamiltonian. Although DFT has been quite successful 
in describing many materials it is not clear a priori, whether it can give a reasonable 
description of the transport properties of nanoscopic conductors. This issue is discussed 
in Ch. El 

The effect of inelastic scattering by electron-phonon interaction on the electrical trans- 
port through Au nanocontacts [68] and atomic Au wires [75] has been studied recently 
with microscopic models on the basis of the Meir-Wingreen formula with different ap- 
proximations for the electron-phonon coupling. Viljas et al. treated the electron-phonon 
interaction in lowest order perturbation theory in combination with a tight-binding model 
of the electronic structure, while Frederikson et al. treated the interaction within the self- 
consistent Born approximation (SCBA) [74j in combination with an ab initio description 
of the electronic structure on the DFT level. 

In principle, it is also possible to introduce decoherence by inelastic scattering into 
the Landauer formalism phcnomcnologically by introducing "conceptual voltage probes" 



37 



2 Quantum theory of electron transport 

which act as phase-breaking scatterers as was first reahzed by Biittiker[63]. However, this 
can only be done in a meaningful way for very simple models systems j76j . To obtain 
meaningful results for realistic models of nanoscopic conductors it is necessary to develop 
a microscopic theory for the effect of inelastic scattering on the transport like the ones 
mentioned above. 
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Density functional theory (DFT) has become a standard method for electronic structure 
calculations in condensed matter physics and quantum chemistry [77j . This chapter ex- 
plains how the transport theory described in Ch. [2]can be combined with DFT to describe 
transport through nanostructures from first principles. Therefore, I will first briefly in- 
troduce the basic principles of DFT in Sec. 13.11 An extensive review of DFT can be 
found e.g. in Ref. |7S]. Finally, in Sec. 13. 2[ I will describe how the DFT based quantum 
transport approach is implemented in the ALACANT package. 

3.1 Density functional theory 

The electrons in a solid or a molecule are described by a many-body wavefunction "^{{ri}) 
of the electron coordinates {ri}, which is a solution of the many-body Schrodinger equa- 
tion: 

n\'^)^E\^). (3.1) 

Ti. is the Hamiltonian of the interacting electron system, containing the kinetic energy T, 
the Coulomb attraction of the atomic nuclei (ion cores) Ve-i and the electron-electron 
interaction Ve-e- 

n^f + Ve-,+Ve^e- (3.2) 

The many-body Schrodinger equation can in principle be solved by expanding the wave- 
function in a basis of Slater determinants, thus converting it to a problem of diagonalizing 
matrices. However, in practice this can only be done for systems with very few electrons 
as the dimension of the Slater basis grows with the number N of electrons as A^!. In- 
stead, before the advent of DFT the standard approach in solid state physics and in 
quantum chemistry has been to approximate the many-body wave function by a single 
Slater determinant that minimizes the total energy. This is the Hartree-Fock approxi- 
mation (HFA) which is the basis for more refined approaches in quantum chemistry, like 
perturbation theory, or configuration interaction [70| . However, the HFA gives a rather 
crude description of the electronic structure of solids and molecules. The main reason is 
that the single-determinant approach neglects electron correlations, but instead describes 
the electrons as independent of each other interacting via an effective mean field. In 
particular, the description of metallic solids within the HFA is very bad, describing them 
often as insulators or semiconductors. The above mentioned systematic improvement 
of HFA by perturbation theory or configuration interaction on the other hand becomes 
computationally too expensive for larger number of electrons just as the exact diagonal- 
ization approach described above, limiting these methods to systems with a relatively 
small number of electrons. 

An alternative approach is density functional theory (DFT) where the basic variable is 
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the electron density instead of the fuU many-body electron wave function '^{{fi}): 

n{r) = {r\P\r) =N j |*(r, fa, . . . fjv)|' dfa . . . dfv, (3.3) 
which is a much simpler quantity to handle than the full many-body wavefunction. 

3.1.1 Hohenberg-Kohn theorem 

The basis for DFT is the Hohenberg-Kohn (HK) theorem [IS] which reduces the fully 
interacting A'^-electron problem to determining the ground state electron density n(f): 

1. The non-degenerate ground state energy of an N electron system, Ef^^ is a unique, 
universal functional of the electron density n{r). 

2. The electron density that minimizes the energy functional is the exact ground state 
electron density 

So if we know the energy functional of the electron density E[n] then we can determine 
the electron density of the ground state by simply minimizing the energy functional, 
which would be much simpler than resolving the many-body Schrodinger equation. It 
is important to note, that no approximations has been made so far. Thus unlike the 
HFA, DFT is in principle an exact theory which allows us to calculate the exact ground 
state density of an interacting electron system. However, the drawback is that the exact 
energy functional is not known. So one has to come up with approximations for the 
energy functional. Different approximations to the energy functional are discussed after 
the next subsection. 

3.1.2 Kohn-Sham equations 

The HK theorem actually does not give a practical recipe for calculating the energy func- 
tional of the electron density, but only states the existence of that functional. A practical 
way for calculating the energy functional is given by the Kohn-Sham (KS) equations [79j 
which map the interacting electron system with some electron density n{r) onto an aux- 
iliary non-interacting system with the same electron density as the interacting system. 
The KS equations are thus the starting point for any practical implementation of DFT. 
We start by decomposing the energy functional into different contributions: 

E[n] = T[n] + Enln] + E,,[n] + E,,t[n]. (3.4) 

T[n\ is the kinetic energy for that electron density 

T[n] = (*|t|*), (3.5) 

where is the many-electron wavefunction corresponding to the electron density n for 
that system. Note, that there is no explicit formula for computing the kinetic energy from 
the electron density n. This is one of the reasons why there is no explicit formula for the 
energy functional. Eh [n] is the Hartree term describing the classical Coulomb repulsion 
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of the electron cloud. For this term an explicit formula can be written for its functional 
dependence on the electron density: 

i!^2;c[n] is the so-called exchange-correlation (XC) term and no explicit exact formula exists 
for this term cither. It is thus the other term besides the kinetic energy term responsible 
for the impossibility of finding an explicit expression for the total energy functional. The 
XC term contains pure quantum effects like the exchange interaction due to the Pauli 
principle and electron correlation effects. This term is usually expressed in terms of an 
integral over the electron density and an unknown XC energy density exc[n]{r): 

Exc[n] = / drn{f)exc[n]{r). (3.7) 



Finally, the last term Eext[n\ is the interaction of the electron cloud with the external 
potential due to the atomic nuclei. This term can again be written as an explicit functional 
of the electron density: 

EexM - j drn{r)Ve^,{r) (3.8) 

The basic idea of the KS method is to introduce an auxiliary set of one-electron wave- 
functions {(/'z(r)} that give rise to the same electron density n{r) as the full many-body 
wavefunction \I'(fi, . . . , rV): 

N 

i=l 

We can then define a KS energy functional TKs[n] that can easily be calculated from 
the auxiliary KS wavefunctions (/j^: 

N 

TKs[n]=Y.{<j>^\f\(|>,). (3.10) 

i=l 

However, the KS kinetic energy term is not identical with the kinetic energy term of 
the interacting electron system T[n\. The difference between the KS kinetic energy and 
the true kinetic energy is again unknown, and is absorbed into the also unknown XC 
functional which is redefined in the KS method as: 

Exc N -> E'xc [n] = Exc N +T[n]~ Tks [n] . (3.11) 

Thus all many-body effects have now been shifted to the XC functional, which is the 
only contribution to the total energy functional that needs to be approximated. All 
other terms can be calculated exactly within the KS method. The main advantage of 
the KS method is that we get a set of effective one-body Schrodinger equations -the KS 
equations- for the auxiliary KS wavefunctions (pi which in turn give the ground state 
electron density of the true interacting electron system. From the variational principle 
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it follows that in the ground state the total energy must be stationary with respect to 
variations of the KS wavefunctions: 

^-e,0.=O (3.12) 

where the are Lagrange multipliers which ensures the orthogonality of the KS wavefunc- 
tions (j>i and will give the effective KS eigenenergies. The KS equations follow directly: 

+ VH[n]ir) + K.xt(r) + Vxc[n](r)) 0,(r) = e.0,(f) (3.13) 



2m 



where 



VH[n]ir) := / dr^ (3.14) 



is the Hartree potential due to the direct electron-electron interaction and 

V.clnm := '-^(r). (3.15) 
on 

is the effective exchange-correlation potential. 

The interacting problem has now been reduced to a set of non-interacting Schrodinger 
equations. One should emphasize that up until now no approximations have been made. 
The KS equations are exact if one knows the exact XC functional. The KS Hamiltonian 

Hks = -^y^ + Vecct + VhIu] + Vxc[n] (3.16) 

depends on the electron density 7i(r) which is the quantity to be determined by the KS 
wavefunctions. Thus the KS equations are a non-linear eigenvalue problem which has 
to be solved self-consistently. The central task of the KS formulation of DFT is to find 
suitable approximations for the unknown XC functional. 



3.1.3 Energy functionals 

Local density approximation (LDA). The simplest approximation is to start from the 
homogeneous electron gas, and assume that the electron density is only slightly modulated 
by the potential of the ion cores. In this case the XC energy density excinlir) at point 
r can be assumed to be that of a homogenous electron gas with the same density as the 
local density of the non-homogeneous electron gas, i.e. exc'^W(^) ^xc''('^(^))- 

Ekc>] nif)e^^^in{r))dr. (3.17) 

Thus in LDA the XC energy density Cxc"^ is simply a function of the local electron 
density n(r) and not a functional. In spite of its apparent simplicity LDA has been quite 
successful in describing metallic systems. One of its major shortcomings is that it cannot 
give a reasonable description of strongly localized electrons like e.g. the d- and /-states 
of transition metals. In this case the electron density is strongly modulated and the 
approximation of a slowly varying density is obviously a bad one. 
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Generalized Gradient Approximation (GGA). A logical step for improving LDA 
is to make the XC energy density also a function of the gradient of the electron density 
in addition to the local electron density: 



While LDA works well for simple metals, GGA improves the description of transition 
metals considerably. However, in general even GGA does not give a good description of 
semiconductors and insulators. And in particular, it fails to give a good description of the 
so-called strongly correlated materials like the transition metal oxides. This is due to the 
insufficient cancellation of the self-interaction error by the approximate XC functionals 
of both LDA and GGA (although the self-interaction error is smaller for GGA). This 
self-interaction error raises artificially the energy of occupied localized states like the d- 
and /-states which are of outmost importance in the transition metal oxides [80j . 

Hybrid functionals. As said before the main problem of the LDA and GGA func- 
tionals is the spurious self-interaction error for localized electron states like e.g. the d- 
and /-states of transition metals, and the orbitals of molecules. On the other hand in 
the HFA the exchange interaction term cancels exactly the self-interaction present in the 
Hartree-term. Thus a possible way of correcting the spurious self-interaction in the GGA 
functional is to reintroduce Hartree-Fock exchange (HEX) into the KS theory by mixing 
some HEX with the GGA exchange functional, an idea that was first put forward by 
BeckejST]. Although a deeper theoretical justification was given for this ad hoc correc- 
tion, the exact amount of HEX that needs to be reintroduced in order to give a sufficient 
correction of the self-interaction error in the GGA functional is unfortunately unknown 
and thus needs to be fitted to experiments |81j. Thus strictly speaking DET using hybrid 
functionals does not really represent an ab-initio method. However, it turns out that 
the parameterization of the B3LYP hybrid functional which was determined by fitting to 
experimental results of a certain class of molecules also works surprisingly well for many 
molecules which had not been included into the fitting data set. Moreover, it also gives 
a reasonable description of the electronic structure and of the magnitude of the band 
gap of insulating solids (where GGA normally fails) and even of some strongly correlated 
materials like e.g. NiO [82] . 

LDA+U method. The idea behind the LDA-I-U method [83] is similar to the hybrid 
functional method: Eor the strongly localized atomic orbitals of a material e.g. the d- 
orbitals of the transition metals and the /-orbitals of the rare earth metals, a Hubbard U 
parameter presenting the screened electron-electron interaction in these orbitals is added 
to the LDA KS Hamiltonian and treated in the HEA thus correcting the self-interaction 
error of these orbitals. In principle it is possible to extract the U parameter from a previous 
LDA calculation. However, in practice this is seldom done since it does not reproduce 
the experiments sufficiently well. Instead the U parameter is normally fitted in order to 
reproduce certain material properties just as the amount of HEX in the hybrid functionals 
has been fitted to reproduce certain properties of molecules. As the hybrid functionals, 
LDA+U gives a reasonable description of the electronic and magnetic structure of strongly 
correlated materials like NiO, and improves considerably the magnitude of the band gap 
of these materials in comparison with LDA [M]. 

Exact Exchange functionals. Another approach to correct the self-interaction error 
of standard DET methods is to derive a local KS exchange potential from the non-local 




(3.18) 
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exact exchange energy of the HFA. The Exact Exchange (EXX) scheme is sclf-intcraction 
free since the exact exchange energy of the HFA cancels exactly the self-interaction of the 
Hartree-term. On the other hand, the EXX potential is local as required from the KS 
theory in contrast to the HEX term added by hand in the Hybrid functional approach 
which is inherently non-local. Therefore the EXX together with a good approximation 
of the correlation energy (as is the case in LDA and GGA) comes very close to an exact 
implementation of DFT within the KS scheme, and thus gives an excellent description of 
the electron density. Moreover, no empirical parameters are needed so that EXX is a real 
ab initio method. Unfortunately, the computational cost for the calculation of the EXX 
functional is very high compared to other methods [53 [HSl [S3 • 

3.2 Kohn-Sham based NEGF formalism 

We now extend the Landauer formalism presented in Ch. [H to the case of interacting 
electrons. The system is described in an effective single-particle picture by treating the 
electron-electron interaction by a mean-field method like the KS scheme or the HFA. 
The problem treated here, is to apply a mean-field method (established only for finite 
or periodic systems) to the open (i.e. infinite and non-periodic) system of the transport 
problem consisting of the leads and the device. The arguments are presented for the KS 
formalism, but are easily applied to the HFA. 

Far away from the scattering region the electronic structure of the leads has relaxed to 
that of the bulk (perfect) lead. So if the device region is sufficiently big (i.e. a sufficiently 
big part of the electrodes is included into the device region) the lead electrons can be 
described as non-interacting quasi-particles moving in the effective potential landscape 
(generated by the ion-core potentials and the effective mean-field potential of the other 
electrons of the bulk leads. Thus we describe the leads by the fixed effective one-body 
Hamiltonians of the bulk leads, and the electron-electron interaction is taken into account 
explicitly only the device region. The (many-body) device Hamiltonian Hu comprises a 
single-body term describing the kinetic energy Td and the external potential Ve-i inside 
the device, and the electron-electron interaction Ve-e between the electrons. Applying 
the KS scheme the resulting KS Hamiltonian for the entire system is given by 

/ Hi Hld \ 

UksM = Hi^N Hi,,^ , (3.19) 

V B.RD Hr J 

and depend on the electron density n of the entire system. In order to keep things simple 
we will not treat overlap between basis functions explicitly. The arguments presented 
here are easily generalized to non-orthogonal basis sets by using the formulas presented 
in Ch. [2 

Wc define a Green's function (GF) matrix corresponding to the KS Hamiltonian which 
thus also becomes a functional of the electron density n{r) 

GKs[n]{E) = {E- nKs[n] + ifj^+r'- (3.20) 

On the other hand, we can obtain the electron density by integrating the GF up to the 
Fermi energy: 

n(f) = --Im / dE{f\GKs[n]{E)\r). (3.21) 

^ J —OO 
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Thus, by self-consistently calculating the density matrix from the Green's function, we 
solve the mean-field problem and thus minimize the energy of the entire system. 

Up until now we have only restated the KS formalism in terms of GFs instead of 
wavefunctions. Now we have to apply it to the situation of the open system described 
by the Hamiltonian p.l9p . The difficulty lies in determining the Fermi energy ep for the 
open system. However, as we will see now, this is only a conceptual problem, that can 
solved in a relatively easy manner. 

The Fermi energy is that of the entire system in equilibrium and in general will be 
different from the Fermi energy of any of the two (isolated) bulk leads due to the electric 
field caused by charges present in the device region. However, this electric field will 
only cause a shift in the electrostatic potential of the leads if the device region is chosen 
sufficiently big, since only the direct Coulomb interaction (i.e. classical electrostatic 
interaction) between charges is long range while the pure quantum mechanical part of 
the Coulomb interaction (i.e. exchange and correlation contributions) is relatively short 
range. On the other hand the Fermi level must be the same throughout the entire system. 
Thus we can shift the Hamiltonian of the two leads to a common arbitrary Fermi level 
which we choose to be zero for convenience. The electrostatic shift of the two leads is 
equivalent to a shift of the Hamiltonian of the device region: 

GoMiE) = {E- HdH -Ad- ^l{E) - ^r{E))-\ (3.22) 

The negative of that shift corresponds to the real Fermi energy of the entire system: 
tp = — Ajj. Next, in order to find the real Fermi level of the system wc impose charge 
neutrality on the device region since we have assumed that the device contains a sufficient 
part of the electrodes so that it has relaxed to its bulk electronic structure. The number 
of electrons in the device region corresponding to a certain energy shift is found by 
integrating the Greens function of the device part up to the adjusted Fermi energy of the 
leads, i.e. zero: 

Nd{Ad) = Im / dETT[GD[n]{E)]. (3.23) 

Thus by imposing charge neutrality on the device region we obtain the energy shift A^ 
of the device Hamiltonian and thus the Fermi energy of the entire system. 

The KS Hamiltonian of the device is given by the single-body term = T^ks + 
Vext and the effective single-particle potential of the electron-electron interaction com- 
prising the Hartree- and the XC term Vc[n] = V/f [n] +V[n] which are both functionals 
of the electron density n{f). But since we have assumed that the electrons do not inter- 
act with each other outside the device region the Hamiltonian and consequently the GF 
only depend on the electron density n]j{r) inside the device region. In turn p-o can be 
calculated by integrating the device part of the GF: 

1 

nD{r) = Im / dE {f\G D[n]{E)\r) . (3.24) 

^ J — oo 

The self-consistent calculation of the electron density from the GF to minimize the energy 
of the whole system can thus be restricted to the device's subspace. 

The energy of the system can be calculated by summing up the KS energies corrected 
for the double-counting: 

E[n] = I] - A^d] + E^,[nD] - j dr 14e[ni5](r) n,5(r). (3.25) 
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where we have made use of the fact that the electrons are assumed to be interacting only 
inside the device region so that only there we need to account for the double-counting 
error made in summing up the KS energies. In order to separate the energy contributions 
of the leads from the energy contribution of the device in the sum we rewrite the sum over 
the KS energies in terms of the KS Hamiltonian and the corresponding density matrix P 
defined by the KS eigenstates: 

E'[Pks] I^e. = Trp/^sH^s]- (3.26) 

i 

By dividing the KS Hamiltonian and the density matrix into sub-matrices as before we 
can separate the different contributions to the energy, and write the total energy as 

E'[Pks] =El + Eld + Ed + Ebd + Er, (3.27) 

where we have defined the energy of the device, 

Ed = Tr[Pc H^] 

the energies of the left lead and right lead. 

El = TrpL H^] and Er[Pj,] = Tr[P,jHfl] 

and the coupling energies between left lead and device. 

Eld = TrpLD Hzji] + Tt[Pdl Hld] - 2Rc (TrpLu Udl]) , (3.30) 

and right lead and device 

Erd - Tr[P RD Hdr] + Tt[Pdr 'Rrd] = 2Re (TviP rd Hdr]) . (3.31) 

Since we have assumed that the device region is sufficiently big so that the electronic 
structure in the leads has relaxed to that of the bulk, P^, and Pr remain constant during 
the self-consistent procedure, and thus also the energy contributions of the leads, E]^ and 
E^. On the other hand the coupling energies actually do change because the GFs Gld 
and Grd directly depend on Gd (see eqs. (|B.6p and (|B.7|) ). So the two energy terms 
£^LD and -Erd always have to be included into the energy calculation. 

The arguments presented in this section are easily adapted to the HFA. Instead of the 
KS Hamiltonian Hks the Fock matrix F is used to calculate the GF of the device which 
depends on the density matrix P which in turn can be calculated from the Hartree-Fock 
GF of the device thus defining again an iterative procedure for the self-consistent solution 
of the HF equation. In the HFA the expression for calculating the total energy from 
summing up the HF cigcnencrgies corrected for the double-counting is much simpler than 
the corresponding expression in the KS scheme. 

Etot = ^ - Itt[PdVhf] = Tr[PF] - ^Tt[PdVhf], (3.32) 

i 

where Yhf is the effective Coulomb interaction in the HFA. 

Finally, Fig. 13. II shows a schematic picture of the implementation of the self-consistent 
transport method presented in this section. This approach is implemented in the ALACANT 
package which interfaces the quantum chemistry program GAUSSIANjlT] to implement 
the NEGF formalism. 



(3.28) 



(3.29) 
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Figure 3.1: Diagram illustrating the ALACANT self-consistent procedure for KS based 
NEGF formalism as explained in the text. 
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3.2.1 Modeling of the bulk electrodes 

As explained in Ch. [21 part of the electrodes has to be included in the device region, and 
in the case of nanocontacts the nanoscopic device actually consists only of the atomically 
sharp tips of the two metal electrodes. In the KS based NEGF described before and 
implemented in the ALACANT package the electronic structure of the device region and 
thus the part of the electrodes included in the device is calculated self-consistently within 
the KS scheme while the electronic structure of the rest of the bulk electrodes is assumed 
to be fixed. Their effect on the electronic structure of the device region is taken into 
account via self-energies in the KS GF of the device region, eq. (|3.22p . 

Since the exact atomic structure of the macroscopic electrodes in a real experiment 
is unknown the question arises how to model the part of the left and right semi-infinite 
electrodes that has not been included in the device region, and how to calculate the 
corresponding self-energies. The ignorance of the exact atomic structure of the electrodes 
in an experiment is a not negligible source of uncertainty in the transport calculations. 
Therefore an appropriate model for the bulk electrodes should have the property that 
the calculated conductance does not depend too much on the actual atomic structure of 
the electrode. A first problem that one encounters is that the self-energy of the bulk 
electrodes can only be calculate for idealized situations. 

A possible choice is to model the electrodes as perfect semi-infinite crystalline leads 
of some finite thickness, i.e. nanowires |SH1 1131 [SH]- In this case the self-energies can 
be calculated as shown in Ch. [5] and App. [C] However, it has been shown that the 
results of a transport calculation using nanowires as electrode models depend strongly on 
the actual thickness of the nanowire and the actual atomic structure [SOI EI]- Moreover, 
the electrodes in real experiments are very far from being perfect nanowires but are 
actually substantial polycrystalline electrodes. Nevertheless, using perfect nanowires as 
the electrodes is interesting for performing model calculations, e.g. to study the effect of 
a single impurity or vacancy on the transport properties of an otherwise perfect nanotube 
or the effect of a constriction in a perfect graphene nanoribbon [52] . 

In order to avoid the that the resulting conductance reflects the finite thickness of 
the nanowire one can go to the crystalline limit by making the wire infinitely thick. In 
practice this can be done by making the system periodic in the directions perpendicular 
to the transport direction so that the self-energy matrices become dependent on the 
wave vector perpendicular to the transport direction |911 154j . However, this approach is 
computationally expensive and the resulting conductance still reflects the crystal direction 
of the electrodes. Alternatively, jellium models have been used in methods based on a 
description in terms of scattering states [Ml [M|- In this case the electrodes are completely 
featureless and thus do not reflect any of the chemical properties of the material. Thus 
this choice is also not free from controversy [55] . 

For these reasons the choice in the original implementation of the ALACANT package 
(then called GECM for Gaussian Embedded Cluster Method) [52] was to describe the bulk 
electrodes with a Bethe lattice parameterized tight-binding model with the coordination 
and parameters appropriate for the chosen electrode material. The advantage of choosing 
a Bethe lattice (BL) resides in that, although it does not have long-range order, the short- 
range order is captured and it reproduces fairly well the bulk density of states of most 
commonly used metallic electrodes. The right hand side of Fig. 13.21 depicts schematically 
a BL of coordination 6. The left hand side of that same figure illustrates schematically 
how the device (depicted here as a nanocontact) is connected to the BL electrodes: For 
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Figure 3.2: Left: Schematic 2D illustration of a nanocontact (big grey circles) with the 
first atoms of the BL (small white circles) attached to the outer planes of the nanocontact. 
Right: Finite section of Bethe lattice (BL) with coordination 6. All atoms of the BL have 
the same coordination as in the corresponding crystalline structure giving rise to short 
range order. But there is no long range order in the BL due to the absence of closed 
loops. 

each atom in the outer planes of the device, a branch of the BL is added in the direction 
of any missing bulk atom (including those missing in the same plane). The directions in 
which tree branches are added are indicated by white small circles which represent the 
first atoms of the branch in that direction. Assuming that the most important structural 
details of the electrode are included in the central cluster, the BLs should have no other 
relevance than that of introducing the most generic bulk electrode for a given metal. In 
App. [D]it is explained how to calculate BL self-energies. 

Fig. 13.31 shows the result of an actual transport calculation for the Al nanocontact as 
shown on the left hand side of that same figure. The left hand side of Fig. 13.31 also shows 
how the surface atoms of the BL electrodes arc connected to the outermost atomic planes 
of the nanocontact. The electronic structure of the device region has been calculated 
on the LDA level of DFT and the CRENBS minimal basis set and effective core pscudo 
potential (ECP) by Christiansen and coworkers [96] has been employed. 

Apart from giving a reasonable generic description of the bulk electrodes the computa- 
tion of the self-energies for the BL model is also computationally very efficient compared 
to computing the self-energies for perfect crystalline systems due to the absence of loops. 
In this respect the BL model for the description of electrodes is advantageous as compared 
to the perfect crystalline models of the electrodes. However, if one is interested in study- 
ing transport in an idealized model system in order to understand fundamental aspects 
it might be important to avoid any kind of scattering that occurs at the electrode-device 
interface, e.g. to study the intrinsic transport properties of nanotubes, graphene nanorib- 
bons or metallic nanowires neglecting the otherwise inevitable scattering at the interface 
between the bulk electrode and the nanotube or nanowire. Therefore a part of this the- 
sis has been dedicated to the implementation of a module for calculating self-energies of 
perfect one-dimensional electrodes into the ALACANT package. 

In Ch. [2]and App. Owe have explained how to calculate the self-energies in the case of 
ideal one-dimensional leads of some finite width taking into account the overlap between 
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Figure 3.3: Left: Atomic structure of an Al nanocontact consisting of two perfect pyra- 
mids along the (OOl)-direction made up of 14 Al atoms (shown in grey) each. The surface 
atoms of the BL connected to the outermost planes of the nanocontact are also shown in 
a gold-like color. Right: Transmission of the Al- nanocontact shown on the left. The elec- 
tronic structure of the device region was calculated on the LDA level of DFT employing 
the CRENBS minimal basis set [SS]. The tight-binding parameters of the BL are fitted 
to reproduce LDA electronic structure calculation of crystalline bulk Al and were taken 
from the handbook by Papaconstantopoulos . 

the atomic orbitals of the atoms making up the one-dimensional lead. We now apply the 
developed theory to calculate the transport properties of an Al nanowire (Fig. 13. 4p based 
on previous DFT electronic structure calculation of the nanowire on the LDA level. We 
use the CRYSTAL03 ab initio program for periodic systems [98j and employ the same 
CRENBS minimal basis set with ECP as before [96]. The nanowire is along the 001 
direction of the bulk crystal and the atoms have approximately bulk distances (nearest 
neighbor distance « 2. 8 A). The primitive unit cell of the nanowire (colored section in 
top panel of Fig. 13. 4p consists of two planes containing 5 and 4 atoms, respectively. By 
defining a unit cell for the transport calculations consisting just of one primitive unit cell 
we can take into account interactions between nearest and next-nearest neighbor planes. 
With a unit cell consisting of two two primitive unit cells we can take into account inter- 
plane interactions up to fourth order. 

Fig. 13.41 shows the DOS and transmission of the nanowire in the case that a) only 
interactions between up to next-nearest neighbor layers (red curves) and b) interactions 
between up to fourth nearest-neighbor layers (green curves) are taken into account, and 
we observe that DOS and transmission change when taking into account interactions 
of higher order. However, beyond fourth order the interactions become negligible, so 
that DOS and transmission do not change further when taking into account interactions 
beyond fourth order. 

In order to compare with transport calculations describing the electrodes with the BL 
model, we calculate now the transport properties of a similar nanocontact as above but 
connected to two semi-infinite Al nanowircs instead of the BL bulk electrodes. As in the 
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Figure 3.4: DOS and transmission for the infinite Al-nanowire shown above. Red curves in 
the DOS and transmission have been calculated taking into account interactions between 
nearest neighbor and next-nearest neighbor planes only while green curves have been 
calculated taking into account inter-plane interactions up to the fourth order. 
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Figure 3.5: Transmission of the Al nanocontact shown above. The electronic structure 
of the nanocontact was calculated using a super-cell which is periodically repeated in one 
direction using the CRYSTAL03 package for periodic systems. 
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case of the perfect nanowire we calculate the electronic structure with the CRYSTAL03 
program for periodic systems by defining a super-cell containing the nanocontact as shown 
at the top of Fig. 13.51 The super-cell has been chosen large enough so that the electronic 
structure of the outermost atoms of the super-cell has relaxed to its bulk value (i.e. of 
the perfect nanowire). The device Hamiltonian is then taken as the (on-site) super-cell 
Hamiltonian of the converged KS calculation. For the calculation of the self-energies of the 
semi-infinite leads the Hamiltonian of the unit cell and the coupling between neighboring 
unit cells are extracted from the electronic structure calculation of the perfect infinite 
nanowire. 

From a comparison of Fig. 13.31 and Fig. 13.51 wc sec that although the two transmission 
curves differ from each other in the detail, some overall features in the two curves are 
actually very similar. For example, for energies of 2.5 eV above the Fermi level we observe 
a plateau of about 3 in the transmission and around the Fermi level the transmission is 
roughly one. On the other hand, in the case of the nanowire electrodes we see that the 
transmission function is much more irregular than in the case of the BL electrodes. The 
strong oscillations of the transmission function in the case of the nanowire electrodes are 
interference effects due to the finite width of the electrodes which also dependent very 
strongly on the actual width of the nanowire, so that the transmission curve changes 
strongly in dependence of the nanowire width. Only in the limit of very thick electrodes 
do these interference effects disappear [SO] • On the other hand the transmission function 
in the case of the BL electrodes and is essentially independent of the details of the details 
of the electrode-device interface as long as the interface is not too close to the nanoscopic 
constriction [55]. 

3.3 Beyond density functional theory 

The methodology presented in this chapter for calculating the transport properties of 
nanoscale conductors by combining DFT based electronic structure calculations and 
the NEGF works quite well for simple nanoscale conductors like metallic nanocontacts 
and nanowires, but has difficulties to cope with more complicated nanoscopic systems, 
like molecular conductors. One of the problems clearly is the spurious electronic self- 
interaction inherent in the standard approximations to the density functional like the 
LDA and GGA as discussed in Subsec. 13.1.31 As explained there, this problem can actu- 
ally be cured by semi-empirical methods like using hybrid functionals [81] or the LDA-I-U 
method [53] that correct the spurious self-interaction, but have the disadvantage of intro- 
ducing empirical parameters which have to be fitted to the material properties and thus 
there predictive power is limited. 

A further problem of the DFT based quantum transport approach is that DFT is a 
ground state theory for the electron density while transport properties also involve the 
excited states of a system. Thus even if we knew the exact density functional (which we 
don't) it is not clear that the DFT based transport approach could actually give a good 
description of the transport properties. For the same reason there is no guarantee that 
an exact DFT would yield the magnitude of the band gap of insulators or the HOMO- 
LUMO gap of molecules. Indeed, it seems to be the case that EXX functionals which are 
already very close to an exact implementation of DFT improve only slightly on band gap 
in comparison with GGA. The hybrid functionals and LDA-I-U yield reasonable values of 
the HOMO-LUMO gap and of the band gaps of insulating solids because they go beyond 
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DFT. 

Finally, while DFT is in principle a many-body theory in the sense that DFT yields the 
exact ground state electron density of the correlated many-body ground state if we knew 
the exact density functional, our KS based transport theory is ultimately an effective one- 
body theory as it makes use of the effective one-body KS orbitals and energies. Therefore 
the KS based transport theory presented here does not capture the true many-body effects 
due to electron correlations in the transport properties of nanoscopic conductors like the 
Kondo effect [rTOlfTUT] . 

There are several possibilities to solve these problems. The first is to use the time- 
dependent density functional theory (TDDFT) |f02j which is a time-dependent extension 
of DFT. TDDFT allows to compute the time-dependent electron density of an interacting 
electron system in a time-dependent external potential. Since a time-dependent external 
field gives rise to electron excitations one can obtain the electron densities and energies of 
the excited many-body states from a TDDFT calculation. Thus TDDFT in combination 
with a good approximation for the density functional (i.e. self-interaction free) like e.g. 
EXX should give a very good description of the energy spectrum of metals, semiconduc- 
tors and insulators as well as of molecules. Accordingly, combining the NEGF transport 
formalism with the KS formulation of TDDFT [1031 1104] should improve systematically 
the description of the transport properties of molecular conductors where the transport 
approach based on ground-state DFT with the standard approximations to the function- 
al normally fails. However, TDDFT is computationally considerably more demanding 
than ground state DFT, and so is the EXX method in comparison with the standard ap- 
proximations to DFT as explained before. Consequently, up until now TDDFT transport 
calculations have only been demonstrated for relatively simple model systems. 

Alternatively, one can try to systematically improve the KS GF by many-body tech- 
niques. As explained at the end of Ch. [21 the Landaucr formula becomes exact in the 
limit of small bias voltages and low temperatures even in the presence of strong electron 
correlations if the exact GF calculated from the many-body ground state is known. On 
the other hand, given a good approximation to the energy functional, DFT only guar- 
antees a good approximation for the electron density but not for the GF of the many 
body ground-state. Thus it is not clear whether the KS GF calculated from the effective 
KS Hamiltonian actually gives a good approximation to the exact GF of the many-body 
ground state. But the KS GF is a good starting point for applying many-body techniques. 

A systematic improvement of LDA is the GW perturbation theory where the electron 
self-energy is approximated by the product of the GF G and the screened Coulomb in- 
teraction W. GW gives an excellent description of many metallic systems and improves 
considerably the description of semiconductors and insulators improving e.g. the magni- 
tude of the band gap. Therefore, various groups have recently made an effort to implement 
the GW method for the description of transport in nanoscopic conductors [105|, I106| . 

However, when dealing with strongly correlated materials like e.g. transition metal 
oxides a perturbative approach like the GW method is no longer appropriate and more 
sophisticated methods are needed like e.g. the Dynamical Mean Field Theory (DMFT) 
[107j . The DMFT is based on the observation that electron correlation effects are most im- 
portant for the strongly localized electrons of the d- and /-shells while they can normally 
be neglected for the delocalized s- and p-electrons. Therefore it is a good approximation 
to treat the electrons interactions only locally exact while the interaction with the rest of 
the system can be described on a mean-field level and thus enters only as an effective bath 
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via a self-energy. Thus the DMFT maps the interacting electron problem of an infinite 
system onto the Anderson impurity model with a single interacting site connected to an 
infinite systems jlOlj . In this way the DMFT approach calculates a locally exact GF that 
neglects the long-range correlations of the electrons which are normally unimportant. By 
extending the interacting region of the DMFT over atoms of the complete infinite system 
(Cluster-DMFT) one can in principle make DMFT as exact as necessary. However, the 
cost in computation increases exponentially with the size of the cluster. 
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4 Spin transport 



We are primarily interested in electron transport through magnetic nanocontacts. There- 
fore we need to generalize the above formalism to systems without spin-degeneracy. In- 
troducing the spin-degree of freedom the Hamiltonian and the Green's function (GF) of 
the device are now further subdivided into the following spin-dependent sub-matrices: 



Hd^( ^fl ^f!) G.(i^).f «f!;^) "y^l ] (4.1) 



The magnetization of the electrodes L and R is assumed to be along the z-dircction 
only. Therfore the spin- mixing parts ti and |t of the self-energy matrices of the electrodes 
and hence the coupling matrices become zero: 

Thus the transmission can be resolved into the different spin-channels: 
T{E) = Tv[Tl{E)Q^j,{E)Tr{E)Qd{E)] 

= Tr[(rL(i?) G^oiE) Tr{E) Gd{E)U„,] 

= ^ TT[r2{E){G-^^-^{E))^T'^^^iE)G'S''HE)] 

= T^\E) +T^^{E) +T^^E) +T^^{E). (4.3) 

where we have defined the spin-resolved transmission probabilities 

j.<Ti<T2(^) ■-Tr[T2{E){G''^'''^{E))^T''j^{E)G''J'''{E)]. (4.4) 

rpcri(72 jg j-j^g probability that an electron enters the device from the left electrode with 
spin CTi and is emitted to the right with spin tT2- The meaning of the spin-channels is 
illustrated in Fig. 14.11 

Obviously, when there arc no spin-mixing terms in the device Hamiltonian (hH = 
= 0) then the cross terms T^^ and T-^T 

must be zero, in real nanocontacts, leaving only 
the most important interactions to capture the essential physics of magnetic materials. 




Figure 4.1: Scattering of spin channels. 
An electron with a certain spin-state a 
coming in from the left can be scattered 
to another spin-state a' of the right lead 
or conserve its spin when passing the de- 
vice. 
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Figure 4.2: Energy dispersion for classical scJ-model for homogeneous magnetization and 
for different limits of the exchange splitting J. For J = (left) the chain is paramagnetic 
and for J > it is magnetic. When J becomes larger than half of the band width 2t a 
gap opens between the two spin bands. 

4.1 Domain wall scattering in the classical sd-mode\ 

What is the effect of a domain wall i.e. a magnetization reversal on a spin-polarized 
current? And how docs this effect depend on the concrete magnetization profile? To 
answer this question we consider a tight binding chain with an exchange splitting J for 
the spin, i.e. the so-called classical sd-model: 

^ = -* E - + 4+1.5-) +E^'-^- (4.5) 

where Mi is the magnetic moment of site i which interacts with the spin Si = {ax^ ct^, cFz)® 
This gives rise to an exchange splitting for the spin at that site: 

Mi - Si^ -~J{sm9i cTx + cos9i cr^) ® (4.6) 

To keep things simple we have assumed that the magnetization lies in the xz-plane, i.e. 
the y-component of the of the magnetization is zero, and that the exchange splitting J 
is the same for all sites. 9i is the angle of the local magnetization with the z-axis. For 
a homogeneous magnetization (e.g. Oi — 0) we get two spin-bands that are split by an 
amount 2 J: /| (fc) = —2t cos{ka) =F J. This situation is depicted in Fig. l4.2l for different 
values of the parameter J. For J = we have of course a the situation of a paramagnetic 
metal, while for J > we have a ferromagnet. For the latter we can further distinguish 
two situations. For J < 2t the two spin-bands overlap. This is the typical situation 
in ferromagnetic metals. On the other hand for J > 2t a gap opens between the two 
spin-bands. And the model either describes a magnetic insulator (at half filling) or a 
half-metallic conductor where the conduction electrons are 100% spin-polarized. 
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J = O.&t 



J = 1.5t 




Figure 4.3: Abrupt domain wall in classical sd-model. Top graph: Sketch of model. Bot- 
tom left: Domain wall transmission (solid red line) for J = 0.5t compared to transmission 
of homogeneous ferromagnetic chain for that parameter. Bottom right: The same as left 
graph but for J = 1.5i. 



4.1.1 Abrupt domain wall 



We now discuss the effect of an abrupt domain wall on the transport properties for the 
ferromagnetic phase of the classical sd-model, i.e when J > 0. As depicted in Fig. 
14.31 the magnetization of each of the leads is homogenous, but aligned opposite to each 
other so that there is a sharp magnetization reversal in the device region. However, the 
magnetization axis is the same for all atoms (we choose z, i.e. 6i £ {0, tt}), so that the 
spin-mixing term cx ax in the Hamiltonian (|4.5p is zero. 

When J > 2t there is no overlap between the two spin bands so that the transmission 
is simply zero for all energies since an electron coming in on one of the leads cannot be 
transmitted through the device to the other lead conserving also its spin. But since there 
is no spin-mixing term in the Hamiltonian the transmission probabilities involving a spin 
flip are zero, i.e. the spin is conserved. Thus in the case of the abrupt domain wall there 
can only be a nonzero transmission when the two spin bands do overlap. Note, that 
even in the case where there is an overlap between the two spin-bands the transmission 
becomes zero for energies where the two bands do not overlap as can be seen from the 
two lower panels of Fig. 14.31 Thus for a half-metallic conductor, (i.e. spin-polarization is 
100% near the Fermi level), the effect of a domain wall is maximal, giving theoretically 
an infinite resistance for an abrupt domain wall. This is what makes half-metals the ideal 
materials for spintronics applications. 
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Figure 4.4: Smooth domain walls. The upper panels show 3 domain walls of different 
lengths. Upper left: Domain wall with 3 atoms. Upper center: Domain wall with 5 atoms. 
Upper right: Domain wall with 10 atoms. The magnetization angle 9i changes linearly 
with the atom position i inside the domain. The lower panels show the corresponding 
transmissions of the domain walls (solid red lines) in comparison with the homogeneous 
ferromagnetic chain (dashed black lines) with same parameter J. J = 0.5t in all 3 cases. 



4.1.2 Smooth domain walls 



Instead of the abrupt magnetization reversal the magnetization could change smoothly, 
as illustrated in the top panels of Fig. 14.41 What would be the effect of this on the 
transmission of the domain wall? Due to the smooth change of magnetization in the 
domain wall, the spin-mixing term cx &x will not be zero anymore and therefore the spin- 
flipping channels will open, i.e. T^-^ and 1 will become flnite. Thus one expects that 
the domain wall resistance decreases. This is exactly what happens, as can be seen from 
the bottom panels of Fig. 14.41 which show the transmissions for different lengths of the 
domain walls when the magnetization angle 6i changes linearly with the position. 

We observe that in general the transmission increases with increasing smoothness (i.e. 
increasing domain length) of the domain. In particular, it becomes different from zero 
at energies where there is no overlap between the two spin-bands. This is due to the 
spin-mixing by the non-coUinear magnetization profile of the domain wall which allows 
electrons to spin-flip when crossing the domain wall. The spin-flip process becomes more 
and more adiabatic when increasing the smoothness of the domain wall thus resulting 
in fewer backscattering for incoming electrons. Therefore the perfect transmission of the 
chain with homogenous magnetization is finally recovered for very smooth domain walls 
as we can see in Fig. 14.41 for the case of the domain wall with a length of 10 atoms. 
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4.2 One-dimensional Hubbard model 



The classical sd-modcl discussed in the previous section is a rather simplistic model for a 
magnetic material. There the spin-polarization of the conduction electrons is due to the 
interaction of the electron spin with a classical magnetic moment localized on an atom 
that gives rise to an exchange splitting J of the electrons. However, in reality the spin- 
polarization of the electrons in a ferromagnetic metal does not arise from the interaction 
with a local magnetic moment but is a direct consequence of the Coulomb interaction 
between the electrons and the Pauli exclusion principle. The Coulomb interaction gives 
rise to the so called exchange term when the two spins of the two electrons are both in the 
same direction. This term is usually negative, i.e. associated with an energy gain for the 
electrons. If this energy gain is greater than the loss in kinetic energy which is associated 
with the two electrons occupying the same spin state then usually the electron spin of 
the material will prefer to align parallel. This is the reason why magnetism in material 
normally originates from the c?-electrons which are strongly localized at the atoms so that 
the hopping to the neighbor atoms is rather small minimizing the loss in kinetic energy. 

The most simple model that captures this physics is the Hubbard model. An intro- 
duction to the Hubbard model can be found for example in Refs. [108j . and |109j . The 
Hubbard model consists of only two spin-orbitals per atom which can split in energy 
because of the on-site Coulomb repulsion U. Atoms are coupled only to their nearest 
neighbors by a hopping t, and Coulomb interactions between atoms are neglected. Thus 
the Hamiltonian of the one-dimensional Hubbard model is given by 



T is the so-called hopping term, which describes the kinetic energy Hamiltonian in the 
tight-binding approximation, and V is the on-site Coulomb interaction. As pointed out 
above the Hubbard model is employed to describe d-band materials, i.e. materials whose 
electronic properties are determined by the d-electrons which are responsible for mag- 
netism. 

In the following we will study the formation of domain walls in one-dimensional ferro- 
magnetic chains with the Hubbard model, and their electronic transport properties. In 
order to mimic the situation in nanocontacts where it is expected that a domain wall 
builds up in the atomic-size constriction [39j , we will restrict the formation of the domain 
wall to a finite region of the chain. Outside that region the electronic and magnetic struc- 
ture is assumed to be that of the homogenous chain. The formation of such a constrained 
domain wall will be calculated self-consistcntly following the Hartree-Fock approximation 
(HFA) as explained in the following section. 

4.2.1 Hartree-Fock approximation - the Stoner model 

The HFA to the Hubbard model is also called the Stoner model [110] . Here we will 
treat the one-dimensional Stoner model in the non-coUinear unrestricted HFA (NC-UHF, 
see App. |E} in order to allow the formation non-coUinear magnetization profiles of the 
domain wall. 

In the NC-UHF we find the following Fock matrix for the Hubbard Hamiltonian (|4.7p : 




(4.7) 




(4.8) 
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The first term in the brackets is the Hartree contribution which describes the pure 
Coulomb repulsion, the second term is the Fock contribution describing the Coulomb 
exchange interaction. The Fock terms for cti = fj2 cancel exactly with the Hartree-terms 
oc p^j"' with CT = CTi = (72 which present an unphysical self-interaction of the electrons, 
correcting thereby the self-interaction error. 

The density matrix can be expressed through the creation and annihilation operators: 

P^^' = (4,c,..)*. (4.9) 

nd 



so the Fock-operator corresponding to the Fock matrix (14. 8p can be expressed in 2 
quantization as: 

T = f +UY,{fl^a)n^a - U Y,{iaC^a)* c\,C,^ ■ (4.10) 
i,a i.{7 

In the next subsection we will study the situation of a homogeneous chain in the Stoner 
model. 

4.2.2 Homogeneous Chain 

In an infinite chain all sites are equivalent. Therefore the Stoner model Hamiltonian 
(|4.10p becomes homogeneous, i.e. the expectation values (nio-) arc the same at each site 
hifiitr) = TT-a-. This allows US to choose the quantization axis for the spin, say z, so the 
spin-mixing terms in the Fock operator vanish: 

:f = f + uJ2w,. (4.11) 

i.cr 

We can diagonalize the Fock matrix of the homogenous chain by introducing Bloch 
states, \k) = ^ Ej exp(»fcj)|j): 

.F = ^F'^(fc)cLcfc., (4.12) 

k,(7 

with 

1 

F-^ik) = exp{tkj){Oa\j'\ja) = e(fc) + f/n^. (4.13) 

i=-i 

€k is the kinetic energy of the Bloch wave |fc) (Hamiltonian T): 

e{k) ^2tcos{k). (4.14) 

The two resulting spin-bands (energy F" (k)) have the dispersion relation e(k) of the non- 
interacting chain but arc split by U{nij — n|) in energy. The corresponding DOS V can 
be written as a sum of the spin-up and the spin-down DOS: 

V{e) = V^{e)+V^{e). (4.15) 

where 

if arccos((e - Un^)/2t) ; if |e - Un^\ < 2 \t\ 
V.{e) ^ { (4.16) 
: elsewhere 
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Figure 4.5: Self-consistently calculated magnetization phase diagram for the one- 
dimensional Stoner model (left) and spin-resolved DOS (right) for the ferromagnetic (FM) 
phase (filling n = 0.25, Coulomb parameter J7 = 6) of the onc-dimcnsional Stoner model. 
The dotted vertical line indicates the Fermi level. 



The spin-dependent occupation numbers n| and n-i are given by integrating the spin- 
resolved DOS up to the Fermi energy fj.: 



n„{^) = / deX'a(e) /(e - /i) = - / 
J-oo I" Jv 

f \^i<Uns~2\t\, 
i arccos ~ Un^)/2t) ; |^ - Un^\ <2\t\, 
1 ■,fi>Uns + 2\t\. 



^ <.!7ns+2|t| ^ 

de — arccos ((e — Ung)/2t) /(e — ^) 

Ung-2\t\ 



(4.17) 



We fix the number of electrons n = + per unit cell (also called filling factor), the 
hopping parameter t and the Coulomb parameter U. The self-consistent solution of the 
problem then yields the chemical potential /.t, the occupation numbers n| and n^, and 
thus the spin-density s = n| — nj^ and the total energy per unit cell. 

The left hand side of Fig. l4.5l shows the phase diagram of the self-consistently calculated 
Stoner model for the magnetization (spin per electron) of the chain in dependence of 
the on-site Coulomb repulsion U and the filling n. Obviously there is a ferromagnetic 
(FM) and a paramagnetic (PM) phase. Since we want to study transport in magnetic 
materials we are only interested in the FM phase here. In the FM phase the Stoner 
model is a so-called half-metal, meaning that the electrons are completely spin-polarized. 
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as illustrated on the right hand side of Fig. 14.51 Only one of the two spin bands is 
partially filled, while the other is always empty. The electrons at the Fermi level and thus 
the conduction electrons are 100% spin-polarized. This has important consequences for 
the spin transport. First, at the Fermi level only one transport channel (the majority spin 
channel) is open, leading to a conductivity of 1 Gq. Second, in the case of a domain, the 
MR will be 100% for an abrupt DW since all traversing electrons will be blocked. The 
latter will be discussed in the next section in more detail. 



4.3 DW formation and scattering 

As said in the beginning of this chapter the DW formation will be restricted to a finite 
region of the chain - the device region (D) . We will make use of the self-consistent proce- 
dure described in in Ch. Ofor computing the electronic structure of an open system in the 
HFA. The semi-infinite regions to the left and right of D, i.e. the electrodes L and R, will 
have bulk electronic structure which does not change during the self-consistent calculation 
of the DW formation. In order to give the correct magnetic boundary conditions for the 
formation of a DW the magnetization of L and R will be aligned antiparallel (AP). 
Following (|2.33p , the GF of the device region is given by 



(4.18) 



where Fd is the Fock matrix of the device in the NC-UHF approximation for the Stoner 
model. The self-energy matrices can be subdivided into the spin-resolved self-energy 



matrices H\ ^ and ^ following (|4.2p which can be calculated from electrodes' self- 



energies following (|2.34|2.35p . 

The left and right lead are described by tight-binding chains with a spin-dependent 
energy shift of e^r = Un^. For this case the Dyson equations (|2.4ip and (|2.42p can be 
solved analytically and one obtains for the spin-dependent retarded self-energy of the left 
(1) and right lead (r): 



^/iE'- 



-4*2 



2 



for E < \\2t\\ 



for ll^ll < ||2t|| 



for E > \\2t\\ 



(4.19) 



Wc follow the procedure described in Sec. 13.21 adapted to the HFA to calculate self- 
consistently the formation of a DW in the AP configuration. The magnetization of an 
atom is defined as the total atomic spin normalized to the number of electrons of the 
atom: 

{S{i)) 



M{i) 



n{i) 



(4.20) 



where the expectation value of the atomic spin is obtained from the density-matrix: 

1 



2Re[pJ/] 
2Im[pl/i 



(4.21) 



P 



TT 
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Figure 4.6: Total energy following cq. (|3.27p corrected for double-counting following eq. 
p.25p in dependence of number of steps in the self-consistent DW formation for a domain 
length of (a) 2 atoms and (b) 4 atoms. Filling n = 0.1 and Coulomb parameter U = 3. 



For a "smooth" DW, i.e. a DW where the magnetization angle changes smoothly with the 
position, the initial guess must have non-coUinear magnetization vectors. In the following 
we will start with an initial guess representing a "linear" DW, i.e. a DW where the 
magnetization angle changes linearly with the position. On the other hand if we want to 
calculate the formation of an abrupt domain we have to start with an initial guess with 
only coUinear magnetization vectors, i.e. the magnetization of each atom of the domain 
must be along the same axis (the z-axis) as the magnetization of the electrodes. 

Fig. 14.61 shows the total energy of two constrained smooth DWs of different length (2 
and 4 atoms) after eqs. (|3.25p and (|3.27p for each step in the self-consistent procedure 
until convergence of the density matrix. Obviously, in the case of a domain length of 
only two atoms the energy is not minimized. This is not surprising if we consider that 
for such a short domain the electronic structure of the electrodes has surely not relaxed 
to the bulk electronic structure as we have assumed. In order to minimize the energy 
we would have to recalculate the electronic structure of the electrodes as well, and not 
let it fixed during the self-consistency as we do here. But here we are interested in the 
formation of constrained domains, since in real nanocontacts the DW formation is also 
constrained by the geometry of the contact [35] . Since we have neglected the geometrical 
part of the problem by considering one-dimensional chains we have to constrain the DW 
artificially. This of course leads to a constrained variational problem, so that not the 
energy is minimized but the free energy defined in a proper way. We note that already for 
a domain of 4 atoms (right panel of Fig. 14.61 ) the constrained variational search minimizes 
the total energy of the system. The assumption of a fixed bulk electronic structure in the 
electrodes becomes better with increasing DW length. 

Nevertheless the resulting energies can be compared to each other and to the energy 
of the bulk chain. Thus we can obtain the cost in energy for the formation of a DW of a 
certain length. Fig. l4.7l shows the DW formation energy per atom as a function of the DW 
length in comparison to the energy of the FM bulk chain. The homogeneous chain has of 
course the lowest energy since the ferromagnetic solution favors the parallel alignment of 
the electron spins to increase the gain in exchange energy. As the DW length is increased 
the energy per atom shrinks because locally the ferromagnetic solution is recovered as the 
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Figure 4.7: Total 
energy following eqs. 
(f?:^ and (p:77|) as 
a function of the DW 
length for the self- 
consistent DW for- 
mation. Filling n = 
0.1 and Coulomb pa- 
rameter C/ = 3. 



DW becomes smoother with increasing domain length. 

The smoothing of the DW with increasing domain length is demonstrated in Fig. 14.81 
The upper panels of Fig. 14.81 show the angle 9{i) of the magnetization vector on the 
i-th atom in the device region (=domain) for a domain length of (a) 6 atoms and (b) 10 
atoms, respectively, while the lower panels show the corresponding magnetization vectors. 
In the first place, one notes that the angle 9 for atoms inside the domain (atoms 1 to 6 in 
(a), and 1 to 10 in (b)) depends almost linear on the position, and that this dependence 
approaches a perfect line with increasing domain length. In the second place, one observes 
that the rotation angle AO = 9{i + l) — 9{i) between the last (first) atom of the left (right) 
electrode and the first (last) atom of the domain is bigger than the rotation angles inside 
the domain. This difference between the rotation angles A6 between electrode and device 
and inside the device becomes also smaller with increasing domain length. Finally, we 
mention that also the atomic charge n{i) and the length of the magnetization vector S{i) 
(not shown) varies slightly with the position, and that the deviation from the bulk 
values of the atomic charge and the magnetization becomes also smaller with increasing 
domain length. In summary, the DW becomes smoother with increasing domain length 
approaching a linear DW where the magnetization angle depends linearly on the position 
and as well the magnitude of the magnetization as the atomic charge are constant and 
equal to that of the FM bulk chain thereby locally recovering the FM state. 

In turn the magnetization profile of the DW affects the electron transport through the 
DW. Fig. 14.9b . shows the total transmission as a function of the energy relative to the 
Fermi energy in the P case and for some domain lengths (2,4,10 atoms) in the AP case. 
One observes that the transmission for the AP case deviates globally from the perfect 
transmission in the P case. The deviation shrinks with increasing length of the domain. 
For a domain length of 10 atoms the perfect P transmission is almost recovered. The 
effect is especially strong near the Fermi energy as the close-up in Fig. 14.9b shows. Near 

^n(i) deviates in both directions by less than 2% from the bulk value n = 0.1 for (a) and (b), while S{i) 
is smaller by 4-10% for (a) and by 2-4% for (b) than the bulk magnetization S = n/2. 
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Figure 4.8: DW profiles for (a) 6 atoms and (b) 10 atoms in the device region. The upper 
panel shows the magnetization vector for each atom in the xz-plane, and the lower panel 
the rotation angle 9{i) of the magnetization vector as a function of the atom position i. 
The leftmost (0) and the rightmost (7 in (a) and 11 in (b), respectively) atom belong 
to the left and right electrode, respectively, which have bulk electronic structure. Filling 
n = 0.1 and Coulomb parameter U = 3. 
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Figure 4.9: (a) Total transmission in dependence of the energy (relative to the Fermi 
energy e^^) of one-dimensional chain in the Stoner model for some domain lengths in 
the AP configuration compared to the perfect transmission in the P configuration, (b) 
Magneto-resistance (using def. MRi = (Tp — Tap)/Tp) as a function of the domain length. 
Filling n = 0.1, Coulomb parameter [/ = 3 for both (a) and (b). 
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the Fermi energy only one spin-band contributes to the conduction in the P case, since 
the ferromagnetic chain is a half-metal as explained in the previous section. Thus for an 
abrupt domain the transmission is completely blocked near the Fermi level: An electron 
passing the DW has to change its spin to travel on in the other electrode since only one 
spin-state is allowed near the Fermi level, and the electrodes are oppositely magnetized 
in the AP case. But an abrupt domain does not allow the electrons to change their spin 
because of the absence of spin-mixing terms in the spin-polarized Hamiltonian. On the 
other hand a non-coUinear magnetization of the domain gives rise to spin-mixing terms 
which can rotate the spin state of an incoming electron. The transmission probability 
becomes higher the smaller the rotation angle of the magnetization vector between neigh- 
boring atoms. Therefore in the case of a non-collinear domain of 2 atoms the transmission 
is nonzero but very small compared to the P transmission. Making the domain larger 
the rotation angles decrease - the domain becomes more smooth - and the transmission 
probability creases until recovering the perfect transmission of 1 in the case of an adia- 
batic DW which locally resembles the ferromagnetic solution. As a consequence the MR 
decreases with increasing length of the domain, as demonstrated in Fig. 14.9b . For an 
abrupt domain the MR is 100% as the DW blocks the transmission completely while for 
a very smooth domain the MR goes to zero: The AP conductance becomes the same as 
the P conductance. 

Finally, we will have a look at the spin-resolved conductance channels. Fig. 14. 101 shows 
the contributions of the individual spin-channels (|4.4p to the total transmission. As 
illustrated in Fig. 14.11 T|| (Tx|) gives the probability of an electron coming from the left 
lead with spin-up (spin-down) and being transmitted to the right electrode conserving its 
spin while (^it) probability of an electron entering with spin-up (spin-down) 

and flipping its spin when being transmitted to the right. The tt-and ||-channel do 
not contribute anything to the transmission near the Fermi energy, in accordance with 
the above reasoning that in the half-metallic limit only one spin-band in the electrodes 
contributes to the conductance, so that in the AP case an electron can not conserve its 
spin when being transmitted to the other electrode. On the other hand a non-collinear 
magnetization of the domain opens the f j- and jf-channels which do not contribute to the 
transmission in the P case. But only the ^-channel gives rise to a non-zero transmission 
near the Fermi level because of the electrons near the Fermi energy being completely spin- 
up polarized in the left and spin-down polarized in the right electrodes. The contribution 
of the ti" and |t"Channels to the total transmission grows with increasing domain length as 
the spin-mixing becomes more and more important. For the same reason the contribution 
of the tt-and ||-channel vanish with increasing domain length. Finally, in the limit of an 
adiabatic DW only the t|- and |t-channels contribute and have the same transmission as 
the tt" and ||-channel in the P case illustrating once again how the ferromagnetic case 
is recovered locally for an adiabatic DW. 

4.4 Summary 

In summary, we have calculated the transmission through domain walls in the classical sd- 
model and the Stoner model taking into account non-collinear magnetization profile. For 
the classical sd-model we have assumed a simple linear dependence of the magnetization 
angle with the position inside the domain like in the work by Tatara et al. [Ill] , and 
obtain similar results. Increasing the length of the domain wall its smoothness increases 
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Figure 4.10: Transmission of individual spin-channels in dependence of the energy relative 
to the Fermi level for AP configuration and domain lengths of 2, 4 and 8 atoms. Filling 
n = 0.1, Coulomb parameter U = 3. 
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giving rise to a general increase of the transmission until for very smooth domain walls 
the perfect transmission of the homogenous ferromagnetic chain is recovered. This can 
be understood by the fact that a non- coUinear magnetization allows electrons to spin-flip 
while crossing the domain wall, so that electrons can be transmitted elastically even at 
energies where there is no overlap between spin-bands. In the limit of a perfectly smooth 
domain wall the spin of an incoming electron is flipped adiabatically while crossing the 
domain thus eliminating the backscattering by the magnetic structure. 

On the other hand, we also calculated the formation of a domain wall in the Stoner 
model. We flnd that also here the local magnetization vector changes more smoothly, 
with increasing domain length, although for a finite domain length the DW is not simply 
"linear", i.e. the magnetization angle does not change linearly with the position inside the 
domain as was assumed for the classical sd-model in l4.1l and by Tatara et al. [111] . More- 
over the magnitude of the magnetization vector is not homogeneous for a finite domain 
length, but is modulated along the domain, and thus is different from the magnetization 
of the ferromagnetic bulk. The deviation from the linear DW with respect to the magne- 
tization angles and the modulation of the magnitude of the magnetization becomes bigger 
the smaller the domain, thus introducing extra-scattering in comparison to the perfect 
linear DW. 
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The observation of huge magneto-resistance (MR) in ferromagnetic nanocontacts [28] 
comparable or even exceeding the celebrated GMR effect has given rise to a strong 
interest in these systems over the last years because of its possible implications in the 
context of nanoscale spintronics applications [HI [H [ ^ |501[5T | [5 ^ [55 1 [M 1 [55 1 [5S 1 [55 1 [57]. 
The observed huge MR was explained theoretically by the following reasoning: It is 
expected that the magnetizations of the two sections of the wire which arc connected 
via the atomic-size contact align antiparallel (AP). so that the overall magnetostatic 
energy is reduced. Then a fairly sharp domain-wall (DW) should form at the atomic-size 
contact, since the cost in exchange energy for the formation of the DW is minimized 
due to the low coordination of the contact atom [35]. It has been further argued that 
the sharp DW should give rise to strong ballistic spin scattering resulting in a large 
extra contribution to the ballistic resistance in the AP configuration [40l [41] . Applying 
an external magnetic field the magnetizations of the two wires will align parallcUy (P) 
erasing the DW at the neck, and thus eliminating the extra contribution to resistance. The 
quantity characterizing the effect is the magneto-resistance MR which is defined as the 
difference between the resistances in the AP configuration Rap and the P configuration Rp 
normalized to either Rap (MRi) or Rp (MR2). The definition of MR as MRi {Rap - 
Rp)/Rap is bounded, i.e. MRi < 100% while the second definition of MR as MR2 = 
(i?AP — Rp) / Rp has no upper limit, i.e. MR2 < 00. In the literature both definitions are 
used. Following the above reasoning high values of MR — named "ballistic" MR (BMR) 
for its origin in ballistic electron scattering — in ferromagnetic metal nanocontacts could 
be expected. 

However, BMR in nanocontacts has been a controversial topic since its first observation: 
While some groups have measured huge values of BMR in Ni nanocontacts, comparable 
to or even exceeding GMR values (typically MR2 « 225%) [Ml [3Ql [32l [33l [35], other 
groups have obtained moderate BMR values [3T1 [35] , do not measure any MR effect [31] , 
or even obtain negative BMR values [57] . The question of huge BMR in ferromagnetic 
nanocontacts and its origin is a very important one because of its expected impact on the 
magneto-electronics industry. 

In this chapter electron transport through atomic-size Ni nanocontacts and the effect of 
the formation of a domain wall in the atomic-size constriction on the transport-properties 
is investigated theoretically using ab initio methods. An atomic-size nanocontact is a 
constriction of just one atom in diameter connecting two sections of a wire. Within the 
atomic-size constriction the electrons are confined to effectively one dimension. In order 
to understand the behavior of the electrons in a material at low dimensions it is therefore 
reasonable to investigate the electronic structure of perfect one-dimensional monatomic 
wires before investigating the more complicated situation of a nanocontact with a realistic 
geometry. The investigation of monatomic nanowires of 3d-transition metals with ab initio 
methods has been performed earlier by Smogunov et al. [1121 1113] . As a precursor to 
the work on nanocontacts we have repeated the ab initio electronic structure calculations 
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of Ni chains, and fully reproduce their results. The insight gained by the calculations 
on perfect one-dimensional monatomic wires will help us to analyze and interprete the 
results for nanocontacts with realistic geometries. 

Ni is a 3d-transition metal that crystallizes in the face-centered cubic (FCC) crystal 
structure. The bulk material is ferromagnetic with an atomic magnetic moment of 0.6^b 
which is considerably less than that of his immediate neighbor to the left in the periodic 
table, Co (1.7^b)- The magnetism of Ni originates from the electrons in the d-bands which 
are spin-split. The s-band is not spin-split and thus the s-electrons do not contribute to 
the magnetism. 

The electronic structure of bulk Ni is quite well described within the local spin density 
approximation (LSDA) [114] of density- functional theory (DFT) [115| . Therefore we will 
at first perform ab initio calculations of Ni wires and nanocontacts at the LSDA level 
(Sec. 15. II and Sec. 15. 2p . Then in Sec. 15.41 we will compare the LSDA results with results 
obtained using the hybrid-functional B3LYP which corrects for the self-interaction error 
inherent in LSDA as discussed in Ch. [31 



In order to illuminate the properties of the electronic structure of Ni at low dimensions we 
first study the electronic structure of infinite monatomic chains of Ni using the quantum 
chemistry program CRYSTAL [55], which allows to calculate the electronic structure of 
periodic systems with standard ab initio methods like the Hartree-Fock approximation 
(HFA) and DFT. However, no evidence of chain formation between Ni nanocontacts 
have been reported to date in contrast to Au or Ft. This is also supported by our ab 
initio simulations of breaking Ni contacts presented in the next section. Nevertheless the 
results for the ideal chains will help us to analyze and interprete the more complex results 
obtained for the more realistic contact geometries of the next section. 

First, we optimize the geometry of the chain, i.e. we search for the lattice spacing 
a of the chain which minimizes the energy. Fig. 15.11 shows the energy in dependence 
of the lattice spacing a between 2 A and 2. 2 A obtained with the LSDA functional and 
the CRENBL basis-set with core pseudopotential. The CRENBL [116j effective core 
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Figure 5.2: Electronic band structure with respect to the Fermi energy (zero hne) of 
monatomic Ni chain at equihbrium lattice spacing a = 2.9A calculated with LSDA func- 
tional and CRENBL (a) and DOLL (b) basis set, respectively. The full lines are majority- 
spin bands, while the dotted lines are minority-spin. 

pseudopotential (ECP) is a small core ECP, describing the 10 innermost electrons only. 
The basis set describes the remaining 18 outer electrons of Ni, and is of very high quality. 
The optimum lattice constant of 2.09A is in agreement with that obtained by Wierzbowska 
et al. |117j . and is significantly higher than the nearest-neighbor distance between bulk 
atoms of 2.5A. 

Fig. 15.21 shows the electronic band-structure of the ferromagnetic (FM) Ni chain at 
equilibrium lattice spacing a = 2.09A calculated with the CRENBL ECP basis-set, and 
an all-electron basis-set (DOLL) optimized for Ni bulk |115j . The results for the two 
basis-sets are very similar below the Fermi level and above near the Fermi energy. Only 
for energies well above the Fermi level deviate the band structures appreciably from each 
other. 

The fact that the band structure of the ECP basis-set CRENBL agrees well with the all- 
electron basis-set DOLL below the Fermi energy indicates that the the core electrons are 
very well approximated by the ECP. The reason for the deviation of the band-structures 
for the two basis-sets at higher energies is that the CRENBL basis-set is far larger than 
the DOLL basis-set and thus provides a better description of the outer electrons. Since 
we are above all interested in transport properties which are determined by the electrons 
near the Fermi energy, the deviation of the ban structures well above the Fermi energy is 
not important here. 
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Figure 5.3: Electronic band structure with respect to the Fermi energy (zero hne) of 
monatomic Ni-chain calculated with LSDA functional and CRENBS+P basis set. The 
left panel shows the majority-spin and the right panel the minority-spin bands. 



In a perfect chain no elastic scattering occurs so that all available transport channels 
at some energy transmit perfectly, i.e. have a transmission of 1. The number of channels 
available is given by the number of bands at some energy, so the zero-bias conductance 
is given by the number of bands crossing the Fermi level. Thus in the case of the perfect 
Ni chain there are one majority-spin channel and six minority-spin channels (two of the 
bands are doubly-degenerate, see below) contributing to the conductance. 

The majority-spin and minority-spin bands are split in energy by up to w leV for 
some bands, and the magnetic moment per atom of the Ni chain is ~1.1/ib for both basis 
sets and thus almost twice as high as the magnetic moment of Ni bulk which reflects 
the low coordination of the chain-atoms compared to the bulk atoms: Due to the lower 
coordination of the chain atoms, the bands become narrower than in the bulk so that the 
cost in kinetic energy to pay for the polarization of the bands becomes smaller. This is 
partially counteracted by the decrease in nearest-neighbor distance with respect to the 
bulk which leads to a broadening of the bands. 

Next, we analyze the orbital nature of the different bands of the Ni chain. Fig. 15.31 
shows the band structure of the Ni-chain separately for the two realizations of the spin 
quantum number. In this case the result was obtained with the CRENBS |116j minimal 
basis set with large core (the 18 inmost electrons) ECP. Because of the smaller basis-set 
compared to CRENBL or DOLL and the larger core ECP results obtained with this basis 
are in general less accurate. However, as Fig. l5.3l shows. the changes in the band structure 



74 



5.1 Monatomic Ni chain 



(aid. 
xy-plane 



'3z 



(b) d,, 
xy-plane 



o 



xz-plane 









• 










-plane 


E 


• 




9 



3lE 





Figure 5.4: Sections of the modulus of the (Cartesian) d-orbitals. The (i3z2_^2-orbital 
(a,d) has rotational symmetry around the z-axis (vertical symmetry axis in (d)) and is 
directed along that axis. All other d-orbitals (b,c,e,f) are shown in their respective planes 
of extension. In the direction perpendicular to that plane the extension of these orbitals 
is very small. 



compared to the results obtained with CRENBL and DOLL arc quite moderate, especially 
near the Fermi level. The number of bands crossing the Fermi level remains the same 
(1 spin-up band, 6 spin-down bands) and the bands conserve their overall aspect. Only 
the spin-splitting of some bands has become stronger now (almost 2c V) their bandwidth 
changes slightly and the magnetic moment per atom now is fiB- On the other hand the 
minimal basis-set facilitates the analysis of the band structure a great deal. 

The orbital nature of the different bands is determined by the rotational symmetry 
of the monatomic chain. The eigenstates of the monatomic chain are therefore simul- 
taneously eigenstates of the z-component (for a chain oriented along the z-axis) of the 
angular momentum operator L^. The s-orbital and the (13^2 _,.2 orbital both have angular 
momentum zero, Iz = 0, so there are two bands with angular momentum zero per spin. 
As shown in Fig. 15.31 the spin-bands with angular momentum zero are the broadest of 
the s and d-bands. One is mostly s-type (thus labeled s-band in Fig. 15. 3p with a small 
contribution oid^z^-r'^ while the other is mostly d3j;2_^2-type (labeled (i322„r2-band) with 
a small contribution of s. The s-bands are only slightly spin-split while the (i3z2_^2-bands 
exhibit quite a big spin-splitting of about leV (for CRENBL, Fig. 15.3b ). The s band is 
by far the broadest as expected, since s-orbitals give rise to highly delocalized electrons. 
The orbitals d^z and dyz can be combined to the rotationally symmetric orbitals with 
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Figure 5.5: Effect of electron scattering in the one-dimensional Ni chain. The length 
of the gap d between the two semi-infinite Ni chains is larger than the lattice spacing 
of a ~ 2.lA of the chains and thus gives rise to scattering. The black curves show the 
spin-resolved transmissions in the case of the perfect chain (no scattering) , the red curves 
in the case of d = 1.2a = 2.5A, and the green curves in the case of d = 1.5a = 3.lA. 



angular momentum ±1. Thus two degenerate bands per spin result from the dxz- and 
dy^-orbitals. These are also spin-split by about leV (for CRENBL, Fig. 15.3b .) and are a 
bit narrower than the d322_^2-bands, since the extension of the d^z- and dj^^-orbitals along 
the z-direction is smaller than that of the o?3^2_r2-orbital (see graphical representation of 
d-orbitals in Fig. 15. 4|) giving rise to a stronger localization of the electrons in these bands. 
Finally, the two degenerate bands per spin with angular momentum ±2 are composed of 
linear combinations of the dxy- and (i2.2_j,2-orbitals which are directed perpendicular to 
the z-axis. This leads to a very small overlap of the orbitals in the z-dircction and con- 
sequently to very narrow bands giving rise to strong electron localization in these bands. 
As the other d-bands they are also spin-split by about leV. 

Finally, we study the effect of electron scattering on the conductance of the Ni chain. 
Therefore we calculate the transmission of two semi-infinite atomic Ni chains separated 
by a gap which is larger than the lattice spacing of a w 2.lA of the two semi-infinite 
chains, and compare to the transmission of the perfect chain, Fig. 15.51 We observe 
that the transmission of the spin-down channel near the Fermi level decreases strongly 
as we increase the length of the gap d. The spin-up transmission on the other hand 
is relatively stable near the Fermi level, so that the spin-polarization decreases with 
increasing length of the gap. This can be understood by the fact that the s-clcctrons are 
much less susceptible to scattering than the d-electrons. Therefore the spin-down channel 
which is composed of s- and d-electrons is much stronger affected by the scattering than 
the spin-up channel which consists solely of a single s-type channel. 

In summary, the d-bands are all spin-split by about IcV while the s-band is only 
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Figure 5.6: Left: Ni nanocontact with par- 
allel (P) magnetization of both electrodes. 
Right: Ni nanocontact with antiparallel 
(AP) magnetization of both electrodes. 



slightly spin-split due to its hybridization with the djz^_j.2-haiid. The spin-splitting of 
the d-orbitals gives rise to the filling of the spin-up d-bands while the spin-down d-bands 
are only partially filled. Thus only the spin-down d-bands and both (spin-degenerate) 
s-bands contribute to the conduction. Because of the spin-up rf-bands being completely 
filled and the spin-down d-bands being partially filled the d-bands behave like a half-metal 
and can be modeled by the Stoner model (see Ch. U]). Consequently a sharp domain wall 
in the monatomic Ni chain would lead to a total blocking of the 5 d-channels, while the 
s-channels which behave like a paramagnetic metal are not affected by the domain-wall 
at all. Thus the total transmission of the one-dimensional chain in the AP configuration 
will be 2, while the transmission in the P configuration is 7. Therefore the MR is given 
by (7 - 2)/7 = 5/7 =71% in def. MRi and 5/2 =250% in def. MR2. 

These results for the perfect Ni chain have been obtained earlier by Smogunov et al. 
\1\2\ I113j by a full ah initio calculation of the magnetization reversal in the chain. This 
result could explain some early experiments [28j , but not recent ones j35j . Furthermore, 
to date, no evidence of chain formation in Ni has been reported. Even so, scattering 
at the electrode-chain contact will always be present, and lead to a reduction in the 
transmissions of the d-channels which are very susceptible to scattering. 

5.2 More realistic contact geometries 

Now that we have analyzed the perfect one-dimensional case we are prepared to study Ni 
nanocontacts with more realistic contact geometries. The fact that the electronic structure 
of the one-dimensional chain near the Fermi level is quite well described employing the 
CRENBS minimal basis-set with large-core ECP justifies that we restrict our calculations 
here to that basisQ The bulk electrodes will be described by a Bethe lattice (BL) model 
as described in Ch. [3]which provides a geometry independent description of the electrodes 
with a Bulk DOS which is smoother than the DOS of the perfect bulk crystal and mimics 
an average over both disorder realizations and the actual electrode crystal orientations 

A reference atomic structure of the contact region has been initially taken like that 
shown in Fig. 15.61 Following Viret et al. |31j , we consider the narrowest region to consist 
of two pyramids facing each other, formed along the (001) direction, and with the two 
tip Ni atoms 2.6 A apart forming a dimer. Bulk atomic distances and perfect crystalline 
order are assumed otherwise. Ah initio simulations of the breaking process as the one 
shown in Fig. 15.81 support this choice. We stress that the section of the nanocontacts 
varies in the direction of the current flow. This is the situation in real nanocontacts and 
differs from perfect 1-dimensional systems, studied in Refs. |1121 1113| . and from bulk 
systems studied by Van Hoof et al. jll8j . In this regard, the geometries proposed by 

^We actually performed also calculations with the DOLL and CRENBL basis-set which confirmed that 
the results do not change essentially. 




77 



5 M nanocontacts 




Figure 5.7: Left: Transmission per spin channel of Ni nanocontact for parallel (P) magne- 
tization of both electrodes as shown on left hand side of Fig. 15.61 Right: Same as before 
but for the Ni nanocontact with antiparallel (AP) magnetization of the two electrodes as 
shown on right hand side of Fig. 15.61 



Bagrets et al. |119| are closer to real nanocontacts, but are not backed up by experiments 
or simulations. 

Figure ()5.7p shows the LSDA conductance as a function of energy for both up and down 
spin channels in two situations: (a) Parallel (P) and (b) antiparallel (AP) bulk magnetic 
arrangements. In both cases the self-consistent solution has been forced to respect the 
high symmetry of the nanocontact. In the AP case the self-consistent magnetization 
reverses abruptly between tip atoms. The resulting magnetic moment for the contact 
atoms is ~ I.O^b in both situations. This value is significantly larger than that obtained 
for bulk or surface atoms (w 0.6/iB and reflects the low coordination of the tip atoms 
forming the contact. 

The eigenchannel analysis restricted to the contact atoms as explained in the next 
section reveals that in the P case the majority channel is, for the most part, composed 
of a single sp orbital channel and conducts perfectly around the Fermi energy (set to 
zero) while the minority channel is composed of three orbital channels (one s- and two 
d-like, which conduct roughly the same), and exhibits a transmission strongly dependent 
on the scattering energy. In the AP case the system is invariant under the combined 
transformations that exchange L with R and t with J,, resulting in identical values for the 
conductance of the two spin channels, which now are composed of a dominant s channel 
and a strongly diminished contribution of the d channels. The conductance ratio for this 
particular case is a; = 2.8/3.65 = 0.77. This yields MRi = 23% and MR2 = 30%, which 
is clearly below large MR claims [25} 155]. 

Since the minority channel conductance evaluated at the LSDA level exhibits a strong 
dependence on the scattering energy, we study now whether or not different geometries 
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Figure 5.8: Ab-initio simulation of stretching of atomic-size Ni contact. The outer planes 
(9+9=18 atoms) are displaced in steps of 0.5A starting from a slightly compressed contact 
as compared to the ideal contact geometry in Fig. 15.71 The inner atoms of the nanocontact 
(10 in total) are allowed to relax to local minimum energy in each step of the stretching 
while the outer planes are kept fixed during the relaxation. 

can change the above results qualitatively. In an attempt to explore other realizations of 
the self-consistent potential compatible with the magnetic boundary conditions and the 
experimental information, we perform ab initio structural relaxations as a function of the 
displacement between outer planes in the core cluster. To do so, we consider a cluster 
like that shown in Fig. 15.71 The inner atoms in the cluster (10 in total) are allowed to 
relax to local minimum energy configurations as we stretch. This results, logically, in 
lower energy solutions and in the loss of symmetry, so that the transmission in the AP 
case now becomes slightly spin-dependent. 

In Fig. I5.9l the transmission at the Fermi energy per spin channel for the P and the AP 
configurations are shown as a function of the stretching up to the break-up point, starting 
from a slightly compressed nanocontact. From this figure we see that the conductance 
of the minority channel for the P configuration changes significantly upon small changes 
while the transmission of the majority channel is quite stable against geometrical changes. 
This behavior reflects the fact that the majority channel is mostly s-type, and therefore 
is quite insensitive to geometrical changes. On the other hand the minority channel has a 
strong contribution from the d-orbitals which are very sensitive to geometrical changes due 
to their high directionality (see Fig. 15. 4p . In the AP configuration the transmission of the 
two spin-channels shows an intermediate behavior between the majority and the minority 
channel of the P configuration upon stretching of the contact. This can be understood 
by considering that the two spin-channels are now composed of s- and d-type orbitals, 
where the contribution from the d-orbitals is weaker than in the minority channel of the P 
case, but not negligible as for the majority channel in the P configuration. Consequently, 
the MR, shown in Fig. 15.91 is small and barely changes as the nanocontact is stretched, 
oscillating between 14% and 23%. 

5.3 Orbital eigenchannel analysis 

It is often useful to decompose the total conductance into the contributions of the trans- 
port eigenchannels, first introduced by Biittiker|63|. These arc defined as the linear 



79 



5 M nanocontacts 



Transmission / spin cliannel IVIR [%] 




Displacement (A) Displacement (A) 



Figure 5.9: LSDA transmission per spin channel at the Fermi energy for P (black) and 
AP (blue) configuration in dependence of the stretching of the contact (Fig. 15. 8p . 



combinations of the incoming modes in a lead that do not mix upon reflection on the 
scattering region and present a unique transmission value [3S]. The decomposition of 
the measurable total transmission in terms of the transmissions of these eigenchannels 
simplifies considerably the interpretation of the results. Knowledge and analysis of the 
eigenchannel wavefunctions would, in turn, allow one to make predictions regarding the 
behavior of the conductance upon distortions of the geometry or other perturbations of 
the scattering region[53j- Unfortunately, in the NEGF approach it is not straightforward 
to extract the orbital composition of the transport eigenchannels. Only the eigenchan- 
nel transmissions can be obtained easily in the NEGF approach from the non-negligible 
eigenvalues of the transmission matrix. However, the associated eigenvectors turn out to 
be useless as obtained. The reason is that these eigenvectors contain the contributions 
to the eigenchannel wavefunctions of the atomic orbitals at one of the borders of the 
scattering region immediately connected to the leads, but not inside. 

In this chapter we present a method for analyzing the orbital contributions to the trans- 
port eigenchannels at an arbitrary cross-section of a nanoscopic conductor by calculating 
the transmission matrix projected onto that cross-section. Our approach generalizes pre- 
vious work by Cuevas et al. |120| for tight-binding-type Hamiltonians to non-orthogonal 
atomic orbitals basis sets as those commonly used in quantum chemistry packages. An 
alternative approach to investigate the contributions of certain atomic or molecular or- 
bitals to the conductance consists in directly removing the respective orbitals from the 
basis set |121| . 

For the sake of clarity we will repeat here some of the important formulas of the NEGF 
presented in Ch. [H Fig. 15.101 shows a sketch of a nano-constriction connecting two bulk 
leads. We assume that the leads are coupled only to the constriction but not to each 
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Figure 5.10: Sketch of the scattering prob- 
lem. L: Left lead. D: device. R: right Lead. 
0: cross-section of interest. 
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other. The Hamiltonian describing this situation is then given by the matrix 



H 



Since many density functional theory (DFT) codes work in non-orthogonal basis sets, we 
also allow explicitly for overlap between atomic-orbitals given by the following overlap 
matrix: 

/ Sl Sld \ 
S = Sdl Sd Sdr . (5.2) 
V Srd Sr J 

The standard approach to calculate the conductance is to calculate the self-energies 
of the leads from the Green's functions (GF) of the isolated leads, i.e., for the left lead 
Sl(£;) = (Hdl - ESnL)KL{E){lii^B - -BSld) where gi^iE) = {ESi^ ~ Hl)"1 is the GF 
of the isolated left lead and analogously for the right lead. From this we can calculate 
the GF of the device: 

Gu{E) = {ESj, - Hd - ^l{E) - ^r{E))-\ (5.3) 

which, in turn, allows us to calculate the (hcrmitian) transmission matrix 

TiE) = r^{Ey/-'Gl,{E)Tj,,{E)Go{E)T^{E)'/\ (5.4) 

where Tl = ^(Sl — Sj^) and Fr = ^(Sr — Sj^). Typically the leads are only connected to 
the left and right borders of the device and are sufficiently far away from the scattering 
region so that they can be described by a bulk electronic structure. From the structure of 
eq. (|5.4p it follows that only the sub-matrix of T representing the subspace of the device 
immediately connected to one of the leads are non-zero. Thus the eigenvectors obtained 
by diagonalizing T only contain the atomic orbital contributions to the eigenchannels at 
one border of the device region but not at the center where the resistance is ultimately 
determined. 

To investigate the orbital nature of the eigenchannels at an arbitrary part of the device 
we can simply calculate the transmission matrix associated to this part. By choosing 
this region to be a cross-section, like that indicated in Fig. 15.101 current conservation 
guarantees that the so-calculated conductance is approximately equal to the conductance 
calculated from the transmission matrix of the whole device. We want to emphasize here 
that this is really only approximately true for a Hamiltonian beyond the tight-binding 
approximation since hoppings between atoms on both sides beyond the selected region 
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Figure 5.11: Sketch of Ni nanocontact consisting of 
two pyramids facing each other along the (001) direc- 
tion with the two tip atoms forming a dimer bridge. 
The device region (grey circles) consists of 28 Ni atoms 
and the left and right electrodes (empty circles) are 
modeled by BLs with appropriate tight-binding pa- 
rameters to reproduce Ni Bulk DOS. 



are neglected. Of course this approximation becomes better the thicker the chosen cross- 
section is. We proceed by further subdividing the device region. The cross-section of 
interest will be referred to as while the regions on either side will be denoted as 1 and 
r, respectively: 
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As mentioned above we will neglect the hoppings (and overlaps) between the left and 
right layers outside the region of interest so we set hir = hri = sir = Sri = 0. With this 
approximation the GF matrix of the cross-section can be written as 

Go{E) = (Eso - ho - S;(i?) - KiE))-\ (5.6) 

The self-energy matrices representing the coupling to the left and right lead, S{(i?) = 
(hoi - £;soi)gi(£;)(hio - Esio) and 'S'riE) = (hor - £;sor)gr(£')(hro - £^Sro), are given by 
the GF of the left layer 1 connected only to the left lead L and the right layer r connected 
only to the right lead R, respectively: 

SiiE) - (Ssi-hi-SL(S))-i (5.7) 
gr(S) = {Es,-h,-'En{E))-\ (5.8) 

The reduced transmission matrix (RTM) with respect to the chosen cross-section is now 
given by 

T'(^) = T[{Ey^^GliE)T[.{E)Go{E)T[{Ey^^ (5.9) 

with r; = im - S;^) and r; = i{T,[. - S;^). Diagonahzing T{E) now yields the 
contribution of the atomic orbitals within the cross-section to the eigenchannels. 

In the following we apply the above described method to analyze the orbital nature of 
the conducting channels of Ni nanocontacts which have recently attracted a lot of interest 
because of their apparently high magneto-resistive properties [IHl 131] • We consider the 
nanocontact to consist of two ideal pyramids facing each other along the (001) direction 
and with the two tip atoms being 2.6 A apart. Bulk atomic distances (2.49 A) and perfect 
crystalline order are assumed for each pyramid. Just as in our previous work on Ni 
nanocontacts [55] we perform ab initio quantum transport calculations for this idealized 
geometry. To this end we use our code ALACANT (ALicante Ab initio Computation 
Applied to Nano Transport). The electronic structure is computed on the LSDA level 
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Figure 5.12: Transmission functions of open transport channels for the Ni nanocontact 
sketched in Fig. KU\ as calculated from the FTM T{E) (solid line) and from the RTM 
T'{E) (dashed lines), (a) shows the only contributing M channel and (b)-(c) the three m 
channels. See text for further discussion. 
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Table 5.1: Eigenvectors of the RTM at the Fermi level for the contact sketched in Fig. 
15.111 Each column gives the weights of the atomic orbitals (AO) given in the left column 
on the tip atom in each eigenchannel shown in Fig. 15.121 
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Figure 5.13: Sketch of Ni nanocontact (27 atoms). As 
in Fig. I5.11l thc contact consists of two pyramids along 
the (001) direction but now both pyramids share the 
same atom at the tip. The device region (grey cir- 
cles) consists of 27 Ni atoms and the left and right 
electrodes (empty circles) are modeled by BLs with 
appropriate tight-binding parameters to reproduce Ni 
bulk DOS. 



of DFT with a minimal basis set and the electrodes are described by means of a semi- 
empirical tight-binding BL model. 

As indicated in Fig. 15.111 we calculate the RTM T'{E) for one of the tip atoms of the 
contact (labeled with 0) and diagonalize it to obtain the eigenchannels and the corre- 
sponding transmissions projected on the tip atom. In Fig. 15.121 we compare the individ- 
ual channel transmissions calculated on the one hand from the full transmission matrix 
(FTM) T{E) and on the other hand from the RTM T'{E). Though the electron hopping 
between regions 1 and r of the contact has been neglected in calculating the RTM the so 
calculated channel transmissions approximate very well those calculated using the FTM 
so that it is very easy to relate the RTM channel transmissions with the FTM channel 
transmission. This shows that the hopping between the regions 1 and r on both sides of 
the tip atom is almost negligible. Only for the one majority (M) channel we see a small 
deviation near the Fermi energy indicating that here 2nd neighbor hopping contributes to 
the transmission of that channel. As the eigenvectors of the RTM (see Table [5TT|) reveal, 
this channel is mainly s-type. Since s-electrons are strongly delocalized there is a small 
but finite contribution from second-neighbor hopping explaining the deviation between 
the FTM and RTM transmission in that channel. The first minority (m) channel is also 
mainly s-type but now it is hybridized with d^^^-r^ and Pz orbitals. The other two m 
channels are degenerate and mainly dxz- and dyz-type strongly hybridized with p^- and 
Py-orbitals, respectively. 

As discussed in our previous work[SS] the five d-type transport channels for the m 
electrons available in the perfect Ni chain [TT2] are easily blocked in a contact with a 
realistic geometry like that in Fig. 15.111 because the d-orbitals are very sensitive to 
geometry. We have referred to this as orbital blocking. It is not so surprising that the 
dx2-y2- and da;j,-channel which are very flat bands just touching the Fermi level in the 
perfect chain are easily blocked in a realistic contact geometry. These bands represent 
strongly localized electrons which are easily scattered in geometries with low symmetry. 
Interestingly, even the d322_j,2-channel, which for the perfect chain is a very broad band 
crossing the Fermi level at half band width, does not contribute to the conduction as 
our eigenchannel analysis shows. This channel is blocked because the (i3j,2_^2-orbital 
lying along the symmetry axis of the contact is not "compatible" with the geometry of 
the two pyramids. On the other hand the dxz- and dy z-channels are both open in that 
geometry because their shape is compatible with the pyramid geometry of the contacts. 
This illustrates how the geometry of a contact can effectively block (or open) channels 
composed of very directional orbitals. Of course, for different geometries we can expect 
different channels to be blocked or opened. 
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(a) majority chan. 1 (b) majority clian. 2,3 (c) minority clian. 1 
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Figure 5.14: Transmission functions of open transport channels for the Ni nanocontact 
sketched in Fig. [Oas calculated from the FTM T{E) (sohd lines) and from the RTM 
T'{E) (dashed lines). The RTM transmissions calculated for the cross-section labeled 
with 0' in Fig. 15.3! are given by the thin solid curves (labeled RTM'). 
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Obviously, the approximation made in the calculation of the RTM becomes worse the 
bigger the hopping between the regions 1 and r is. For example, in the contact geometry 
shown in Fig. 15.31 electron hopping from the layers immediately connected to the central 
atom (labeled 0) is certainly bigger than in the geometry of Fig. 15.111 Indeed, Fig. 15.141 
shows that for almost all channels the RTM transmissions differ appreciably from FTM 
transmissions, making it difficult in some cases to relate them to each other. Fortunately, 
we can judge by exclusion which RTM transmission relates to which FTM transmission 
since for the other channels at least the RTM transmission function mimics the overall 
behavior of the FTM transmission function. However, for more complicated situations it 
might be impossible to match the RTM transmission with the FTM transmission for all 
channels. The cure to this problem is obvious: One has to choose a bigger cross-section, 
i.e., add an atomic layer to the cross-section so that the hopping between 1 and r becomes 
small again. If we choose, e.g., the cross-section labeled with 0' in Fig. 15.31 (including the 
atomic layer to the right of the central atom) the so calculated RTM transmissions now 
approximate very well the FTM transmissions as can be seen in Fig. 15.141 

In summary, we have shown how to obtain the orbital contributions to the eigenchan- 
nels at an arbitrary cross-section of a nanoscopic conductor. The method has been imple- 
mented into our ab initio quantum transport program ALACANT and we have illustrated 
the method by exploring the orbital nature of the eigenchannels of a Ni nanocontact. The 
method works very well when the chosen cross-section is thick enough so that hopping 
from the layers left and right to the cross-section becomes negligible. Hence in some 
cases an additional atomic layer has to be included to the cross-section we are actually 
interested in. Taking this into account the method has no limitations and can be readily 
applied to ab initio transport calculations in all types of nanocontacts [55] and molecular 
junctions. 

5.4 The self-interaction problem 

LDA provides a commonly accepted description of the electronic structure of bulk and 
surface ferromagnetism in transition metals. However, the low coordination of the atoms 
in nanocontacts might give rise to a further localization of the d-electrons (compared with 
bulk). But as discussed in Ch. O LDA fails to describe localized electrons properly due 
to the spurious self-interaction. And although GGA improves somewhat on LDA with 
regard to the self-interaction problem by taking into account derivatives of the electron 
density, this is still not sufficient for materials with strongly localized d-electrons like e.g. 
the transition metal oxides. As explained at the end of Sec. I3.1l an alternative approach to 
the electronic structure comes from the use of a hybrid functionals which reintroduce some 
Hartree-Fock exchange (HFX) in order to correct the spurious self-interaction. B3LYP 
for example combines HFX with a certain GGA exchange functional [81], and happens 
to give a reasonable description of the electronic structure and local magnetic moments 
in NiO [52] and La2Cu04 |122] . Popular alternatives to using hybrid functionals are e.g. 
LDA+V [Sg and SIC [mj . 

Here, we explore how the use of different functionals which improve on LDA with 
regard to the spurious self-interaction affects the transport properties of Ni nanocontacts. 
We start by considering again the perfect monatomic Ni chain. Fig. 15.151 shows the 
transmission functions of both spin-channels for different density functionals. We note 
that the transmission of the perfect chain barely changes when using a GGA functional 
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Figure 5.15: Spin resolved transmission function of perfect Ni chain for different func- 
tional: (a) LSDA, (b) GGA, (c) hybrid functional with 10% of HFX and (d) hybrid 
functional with 20% of HFX. Solid lines are majority-spin and dashed lines minority-spin 
transmissions. 



(b) instead of the LSDA functional (a). Only when introducing HFX by using the B3LYP 
hybrid functional does the transmission change appreciably. With B3LYP. the number of 
conduction channels is reduced to 1 majority -I- 4 minority spin-channels at the Fermi level 
compared to the LSDA results. This agrees with recent LSDA-I-U calculations reported 
by Wierzbowska et al. |117j where the two degenerate flat minority bands dxy, d,^2_y2 are 
shifted downwards in energy because of the exchange interaction canceling part of the 
self- interaction of the strongly localized electrons in these flat bands. 

As can be from Fig. I5.16f a') computing the conductance of the Ni nanocontact shown 
in Fig. 15.61 with a GGA functional, again the result does not change very much with 
respect to the LSDA result (compare with Fig. 15. 6p . However, we should note that when 
using GGA more elaborate basis sets than the minimal basis set employed in the previous 
calculations on the LSDA level are required in order to recover the LSDA result. If the 
same minimal basis set as in the LSDA calculations is used with the GGA functional the 
conductance of the d-type minority-channel is considerably reduced. This seems to have 
to do with the fact that due to its dependence on the gradient of the electron density the 
GGA functional is more sensitive to the geometry than the LSDA functional and more 
complete basis sets are required to get reliable results. With B3LYP the results for the 
conductance (see Fig. 15.16( b)) are remarkably different in regard to the minority channel 
conductance which is strongly reduced at the Fermi level. The minority d-channels are 
most strongly affected by this reduction of the transmission. Why the minority-spin 
channel transmission is so strongly suppressed when using the B3LYP functional in the 
case of the nanocontact but not in the case of the perfect infinite chain might have 
to do with the d-electrons of the tip atoms becoming very localized so that they are 
more strongly affected by the HFX, but further investigation of this issue is needed. 
However, by comparison with the Ni conductance histogram [34] one can see that the 
results obtained with LSDA and GGA are in very good agreement with the first peak 
in the histogram at - 2.4 - S.Oe^/Zi while the B3LYP value is not. So LSDA and GGA 
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Figure 5.16: (a) Transmission per spin channel calculated with GGA density functional 
in the P configuration for the model nanocontact shown in Fig. 15.61 (b) same as (a) but 
calculated with B3LYP hybrid functional. 



functionals might be more appropriate for describing purely metallic nanocontacts than 
B3LYP. However, this might be different when oxygen adsorbates are present in the 
contact region, as will be discussed in the next Chapter. Clearly, LSDA and GGA are not 
appropriate for describing the conduction through a molecule contacted by atomic tips 
as they do not reproduce well the energy levels of molecules. On the other hand B3LYP 
does reproduce molecular energy levels extremely well. So the question of the appropriate 
functional to use in DFT based transport calculations is an open one, and it seems that 
depending on the application one has to employ different functionals and justify their use 
by comparison with experiments. This however limits severely the predictive power of 
DFT based transport calculations. A systematic improvement of the DFT based transport 
calculations can only be achieved by employing many-body techniques on top of the DFT 
calculations like the GW perturbation theory or DMFT as explained in Sec. 13.31 



5.5 Discussion of results 

As a rule of thumb, the maximum number of conduction channels in atomic-size con- 
tacts is roughly determined by the number of valence electrons of the contact atom(s) 
|124j . However, as shown above, this hypothetical upper limit is never reached in Ni 
nanocontacts since the transmission of the minority-spin channel being composed of s- 
and d-orbitals is strongly reduced by scattering due the contact geometry. This result is 
impossible to predict without a full atomistic self-consistent calculation. The majority 
channel is s-type. Thus, this channel transmits almost perfectly and evolves smoothly 
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with the stretching of the contact giving a stable contribution of T|| w 1 (sec Fig. I5.9p . 
The s-orbitals in the minority channel are strongly hybridized with c?-orbitals and, there- 
fore, are more sensitive to the contact geometry. The contribution to the conductance of 
the latter, which form narrower bands, disappears with the stretching and disorder, as 
expected (see Fig. 15. 9p . On the other hand, in the AP configuration mostly one s-orbital 
channel per spin contributes. In the AP case the conductance per spin channel lies thus 
in the vicinity of e'^/h, giving « 2e^/h in total and is fairly stable during the last stage 
of the breaking of the nanocontact. 

Comparing the nature of the eigenchannels contributing to the transmission of the 
nanocontact in the P case with the orbital nature of the bands of the onc-dimcnsional 
chain helps to understand why 3 of the c?-channcls arc blocked in the contact which 
transmit perfectly in the onc-dimcnsional case. An cigcnchanncl analysis restricted to 
the tip atom (see Sec. 15. 3p for the symmetric cluster (Fig. 15. 7p reveals that in the P 
case the only majority channel contributing to the conductance is s-type, while of the 
3 minority channels contributing appreciably to the conductance, one is s-type (with a 
small contribution of the (i322_^2-orbital), one is dxz and one is c^y^-type. The three other 
minority channels that contribute in the perfect chain, namely the doubly-degenerate 
fia;y,d^2_j^2-channel and the (i322_^2-channel do not contribute to the conductance. It is 
worth noting that these channels have not "vanished" at the Fermi level, as can be seen 
from the local DOS (not shown), only their transmission has become practically zero. 
The reason for this lies in the geometry of the d-orbitals. As is already clear from the 
discussion of the perfect chain, the dxy- and (i2.2_y2-orbitals have very small overlap in the 
direction of the z-axis which is the main axis of the contact and of the chain. Therefore the 
electrons in this channel are very localized (flat band in the perfect chain) and thus very 
easily scattered by the geometry. On the other hand, it is a bit surprising at flrst glance, 
that the c?322_r2-channcl which is the second-broadest band (and thus the electrons in that 
band are almost as delocalized as the s-clcctrons) in the perfect chain does not contribute 
either to the conductance in the nanocontact. However, in the pyramid-geometry this 
orbital which is directed along the z-direction (Fig. I5.4b ..d) does not connect well to the 
orbitals in the pyramid and thus is geometrically blocked for transmission. On the other 
hand the minority bands of d^z- and dyz-type contribute appreciably to the conductance 
of the nanocontact since these orbitals connect well to the corresponding orbitals in the 
pyramids. Only by distorting the nanocontact upon stretching the geometrical matching 
is destroyed quite rapidly, leading to a rapid decrease in the conductance of the minority 
channel. The reason behind the very small MR values is thus the orbital (or geometric) 
blocking of most of the a priori available minority channels in the P configuration due to 
the non- ideal geometry of the nanocontacts. 

The results presented in this chapter have been confirmed by recent experiments on Ni 
nanocontacts under very controlled conditions [37l[38], and also by recent theoretical work 
employing slightly different methodologies |1251 [54 lll26j . Moreover, the total conductance 
of « 2.5 — 3.5e^/h for the single-atom nanocontact agrees fairly well with the broad peak 
in the Ni conductance histogram [31]. So what else could be the reason behind the 
huge MR values observed in some experiments when ballistic domain wall scattering 
can be excluded? First, as will be shown in the next chapter, the BMR effect can be 
enhanced considerably (by more than one order of magnitude) due to the presence of 
oxygen surface adsorbates in the contact region [SS]. Indeed, just the experiments that 
obtain huge MR values have not been performed in ultra-high vacuum conditions in 
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contrast to the more recent ones that do not obtain high MR [34l [3ll [38] . But this cannot 
explain the highest MR results obtained [521 [331 [35] which are several orders of magnitude 
higher. Another likely explanation for very high MR values are magnetostriction and 
magnetostatic effects which lead to a distortion of the contact geometry in dependence of 
the relative magnetization of the two sections [3S1 1127j . 
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In going from bulk to lower dimensions material properties often change drastically. A 
recent example is that of interfaces between different insulators which can become metallic 
jl28|, I129j . Even more recently, it has been predicted theoretically that certain oxygen 
surfaces of some insulating ceramic oxides can exhibit magnetism and half-metallicity 
|130j . The ultimate limit in this respect can be found in atomic chains formed in metallic 
nanocontacts which allow to study the transport properties of one-dimensional systems 
of atomic size |131| . Due to the lower coordination of the atoms the properties of metallic 
atomic chains formed in nanocontacts can be remarkably different from those in the bulk. 

However, not all metals form atomic chains in nanocontacts, although recently, it has 
been found that the presence of oxygen favours their formation |132| . For example, 
experiments with Ni nanocontacts [151 [3S1 |3TJ [Ml [351 [3Z] have never shown evidence of 
chain formation. Nevertheless, the presence of oxygen in the contact region could possibly 
lead to the formation of NiO chains. In this context it has also been proposed that the 
rather moderate magnetoresistive properties of pure Ni nanocontacts [55l 11331 11251 154] 
could be enhanced considerably by the presence of oxygen adsorbates on the surface of the 
Ni electrodes |134| . On the other hand, bulk NiO is a common example of a correlated 
insulator with antiferromagnetic (AF) order (see e.g., Ref. |1351 15^). which remains 
insulating even above the Neel temperature when the AF order is lost. Thus it is not at 
all obvious whether or not oxidized Ni nanocontacts or NiO chains should be conductors. 

In this chapter we investigate the electronic and magnetic structure and the transport 
properties of one-dimensional NiO chains, both idealized infinite ones and more realistic 
short ones suspended between Ni nanocontacts. Anticipating our most important results 
our ab initio quantum transport calculations show that short NiO chains suspended be- 
tween Ni nanocontacts can become half- metallic conductors, i.e., carry an almost 100% 
spin-polarized current. This result holds true even for a single O atom in between Ni 
electrodes. Consequently, for antiparallel alignment of the electrode magnetizations the 
transport through the contact is strongly suppressed resulting in very large MR [differ- 
ence in resistance between antiparallel and parallel alignment of the magnetizations of 
the electrodes normalized either to the higher resistance value (MRi) or to the lower one 
(MR2)]: MRi « 90% and MR2 « 700%, respectively. 

6.1 Electronic and magnetic properties of bulk NiO 

Bulk nickel monoxide (NiO) is an example of a strongly correlated material with insulating 
character and antiferromagnetic (AF) order, which has been studied extensively in the 
past, both experimentally and theoretically (see e.g. [HSl [133 [EH [501 [HI [HHl [ESI 
11401 1141] [5^ and Refs. therein). It is now understood that the measured optical gap of 
~4eV is due to charge-transfer excitations from the O 2p band (which overlaps with the 
filled Ni 3d band) to the unfilled Ni 3d band [135| . Thus it is not a pure Mott-Hubbard 
insulator where the gap is solely determined by excitations from filled to unfilled Ni 3d 
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(a) 0% HFX (GGA) (b) 20% HFX (B3LYP) (c) 30% HFX 
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Figure 6.1: DOS of bulk NiO projected onto O atomic orbitals (red lines) and Ni atomic 
orbitals (green lines) for different density functionals. (a) GGA functional, (b) B3LYP 
hybrid functional with 20% HFX, and (c) hybrid functional with 30% HFX. The black 
arrow marks the top of the valence band in each panel. For the calculations the CRYSTAL 
ab initio package has been employed [55] • See text and Ref. [52] for further discussions. 

bands. Though the AF coupling in bulk NiO is strong, it is found experimentally that 
the magnetic order seems not to affect the electronic structure significantly. Even above 
the Neel temperature when the magnetic order is lost, bulk NiO remains an insulator and 
preserves the charge-transfer character and the magnitude of the gap |140j . 

Due to the strong electron correlations of the Ni Bd-electrons, NiO presents a challenge 
for ab initio electronic structure calculations. In fact, the standard approximations of 
DFT — LDA and GGA — do not give a satisfying description of the electronic structure 
of bulk NiO. Because of the insufficient cancellation of the self-interaction in the local 
exchange functional of the LDA, the occupied narrow 3(i-bands are raised in energy. As 
a result LDA sevcrerly underestimates the gap of bulk NiO [5D]. The GGA exchange 
functional improves somewhat the description of NiO but the energy gap is still too small 
and also the charge-transfer character is not captured correctly [80] , as can be seen from 
Fig. EUa). 

As explained in Ch. |31 the self-interaction error inherent in LDA and GGA can be 
corrected by using hybrid functionals like e.g. the B3LYP functional [81] which happens 
to give a reasonable description of the electronic and magnetic structure of bulk Ni [82] . 
As we can see from Fig. IG.lf b). the B3LYP functional which mixes 20% of HFX to 
the GGA exchange functional gives just the right magnitude of the gap of ~4eV. Also 
the charge transfer character of the band gap is somewhat improved although there is 
still a strong contribution of Ni 3c?-states at the upper band edge of the valence band. 
Increasing the amount of HFX beyond the 20% of the B3LYP functional the O 2p-states 
become yet more dominant but also the band gap increases to about 6eV which is above 
the experimental value (Fig. I6.1f c)). With regard to the magnetic structure of NiO, 
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Figure 6.2: Schematic one-electron en- 
ergies of a one-dimensional NiO chain 
in z-direction for minority-spin. To the 
left and right the orbital energies of an 
individual Ni+ cation and anion in 
the crystal field of a one-dimensional 
Ni+O" chain are shown. In the cen- 
ter the formation of valence and con- 
duction bands by hybridization of Ni 3d 
and O 2p orbitals is shown. 

B3LYP predicts the correct antiferromagnetic (AF) order, but slightly underestimates 
the magnetic moments and coupling constants. On the other hand a hybrid functional 
with 30% of HFX gives the correct magnetic moment and coupling constants and also 
predicts the correct AF order. A more detailed discussion can be found in Ref. [52]. Other 
approaches for correcting the self-interaction like the LDA-I-U method the self-interaction- 
correctcd LDA (SIC-LDA), and the GW approximation lead to similiar results as those 
obtained from the hybrid functional approach [83l 11391 [82] . 

For the above electronic structure calculations of bulk NiO we have employed elaborate 
all-electron basis sets for Ni and O |1151 1142] similiar to those employed for reported HF 
and B3LYP calculations of bulk NiO [52], but extended with a diffusive sp-function in 
the case of Ni and a d-polarization function in the case of O which makes them suitable 
also for the description of metallic systems in contrast to the original basis set which have 
been developed for the description of insulating NiO. As can be seen from Fig. 16.11 we 
reproduce previous GGA and B3LYP results for bulk NiO |82] with these basis sets. 

Given the uncertainty with respect to the functional we will employ for the calculation 
of NiO chains in the following both a GGA functional and the B3LYP functional. We 
would like to emphasize here that it is not clear from the outset how much HFX is 
necessary in order to give a reliable description of atomic NiO chains since the electronic 
structure of the one-dimensional NiO chain is quite different from the one of bulk NiO as 
we will see in the following. 

6.2 One-dimensional NiO 

The electronic properties of bulk NiO are, to a large extend, determined by the atomic 
scale properties, like the crystal field splitting of the Ni and O energy levels and the 
amount of electron charge transfered from Ni to O. The latter, in turn, is determined 
by the interplay between Madelung binding energy, the ionization potential of Ni, and 
the electron affinity of O. Due to the lower coordination and the corresponding decrease 
in Madelung binding energy the electron transfer from Ni to O is less favourable in an 
atomic chain than in bulk (where the electron transfer is almost complete resulting in an 
ionic configuration of Ni^+O^"). The proper starting point to discuss the formation of 
energy bands in the one-dimensional NiO chain are therefore the univalent ions Ni+ and 
0". 

In order to understand how the low coordination affects the atomic properties of the 
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Figure 6.3: Energy per 
unit cell of infinite NiO 
chain in dependence of lat- 
tice spacing a calculated 
with B3LYP hybrid func- 
tional. Blue triangles 
indicate the half-metallic 
state (HM), red triangles 
the molecule-like insulat- 
ing state with FM order 
(ML), and black boxes the 
insulating state with AF 
order (AF). 
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constituting ions we have first performed B3LYP calculations of both a single Ni+ ion 
and a single 0~ ion each in the field of point charges that mimic the crystal field of a one- 
dimensional chain of univalent Ni and O ions. We find that for both ions the spin-doublet 
state (5 = 1/2) minimizes the energy as is the case for the free ions. In Fig. 16.21 we 
show schematically the energy levels of Ni+ and 0~ in the presence of the point charges 
for minority spin only. Interestingly, for minority spins, the occupied Ni 3(ii orbitals 
{dxz and dyz) of the Ni+ ion fall energetically in between the occupied and unoccupied 
O 2pi orbital (px and py) of the 0~ ion. Thus the 3di and 2pi orbitals can form two 
filled degenerate bonding bands and two degenerate partially filled antibonding bands as 
indicated in the middle part of Fig. 16.21 The 3d2 {d^y and d^2_yi) doublet is somewhat 
above in energy to the idi but cannot hybridize with the oxygen 2p orbitals. The Ni 3(io 
(d3z2_^2) and 4s are empty while the O 2po {j>z) and 2s orbitals are filled and much lower 
in energy so that no hybridization takes place though symmetry would allow for it. 

On the other hand, for majority-spin electrons (not shown) all five Ni Zd orbitals are 
filled while the 4s is also empty, i.e., the Ni+ valence configuration is ?>d^ and not 4s^3(i^ 
as for the free Ni+ ion. Moreover, all of the O 2p orbitals are filled so that the Ni+ and 
O^ ions can only form either completely filled or completely empty bands for the majority 
spin. Thus the ionic picture suggests that a one-dimensional NiO chain should become a 
half-metallic conductor where only the minority-spin levels form conducting bands. 

Not surprisingly, our calculations for inifinite one-dimensional NiO chains (Fig. 16. 2p 
show that a univalent ionic configuration as an initial guess results in a half-metallic state 
for large separation of the individual chain atoms (i.e., large lattice spacing of 5A) as 
suggested by the ionic picture. However, as can be seen from Fig l6.2r a). the half- metallic 
state is only a metastable state for most values of the lattice spacing. This half-metallic 
state is "shadowed" by a second state with FM order and insulating character. By 
successively decreasing the lattice spacing a of the chain and using the (half- metallic) 
state of the previous step for the initial guess, the half-metallic state can be generated 
also for smaller inter-atomic distances which points towards its metastability. Around the 
equilibrium lattice spacing (a ~ 3.4A) the ground state of the chain has AF order and is 
of insulating character with a substantial gap of ^ 4eV like in bulk. When stretched out 
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Figure 6.4: (a) B3LYP band structure of HM state for a lattice spacing of 3.6A. Solid lines 
indicate majority-spin bands and dashed lines indicate minority-spin bands, (b) Same as 
(a) but for ML FM state. 



of equilibrium the FM state and the AF state become comparable in energy until finally, 
at a lattice spacing of ~ 4.2A, the FM state becomes the ground state. 

The band structure diagram in Fig. I6.4f c) shows that the metastable state with FM 
order corresponds indeed to the half-metallic state suggested by the ionic picture: The 
half-filled doubly-degenerate conduction band is formed by minority-spin Ni 3di orbitals 
hybridized with O 2pi orbitals while the Ni 3d2 orbitals do not hybridize with O 2p 
orbitals and thus form a flat valence band. The lowest-lying empty band is formed by the 
minority-spin Ni 3do orbital which is slightly hybridized with the Ni 4s orbital. On the 
other hand the stable state with FM order and insulating character (sec band structure in 
Fig. I6.4f b)) actually corresponds to the ground state of the NiO molecule which is a ^S"^ 
state |145] . The main difference with the half- metallic state is that now the non- degenerate 
Ni Sdo and 4s orbitals form a minority-spin valence band while the minority-spin doubly 
degenerate half-filled antibonding band composed of Ni 3di and O 2pi bands is emptied 
and a substantial gap of ^ 3eV opens. Thus the infinite chain behaves like an insulator 
for reasonable values of the chain stretching. 

On th other hand using a GGA functional, the one-dimensional NiO chain is always 
conducting for FM order in contrast to the previous B3LYP results. In Fig. 16.61 we show 
the GGA band structure of an ideal one-dimensional NiO chain in the ferromagnetic (FM) 
phase. Compared to the B3LYP band structure for the half-metallic state at same lattic 
spacing shown in Fig. I6.4f a) we see that the occupied bands have been raised consider- 
ably in energy. In particular, the doubly-degenerate flat minority-spin band of type (^2) 
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composed of Ni 3dxy and 3dx2-y2 orbitals (well below the Fermi level with B3LYP), now 
actually crosses the Fermi level. Also the doubly-degenerate and previously half-filled 
minority-spin band of type (di ) composed of Ni id^z and 3dyz orbitals hybridized with 
O 2px and 2py orbital (the only conduction band with B3LYP) calculation are raised 
somewhat in energy. Consequently, the previously empty minority-spin band of type (do) 
composed of Ni Sd^z^-r^ orbitals is lowered in energy with respect to the other filled or 
partially filled 3d bands, and becomes a conduction band. 

On the other hand also the majority-spin bands are raised in energy. The doubly- 
degenerate majority-spin band composed of Ni 3dxz and 3dyz orbitals hybridized with O 
2px and 2py orbitals which was well below the Fermi level in the B3LYP calculation now 
also crosses the Fermi level near to its upper band edge. Thus in GGA the ideal case 
of the infinite NiO chain in the FM phase does not represent a half-metallic conductor, 
although the spin-polarization of the conduction bands is quite strong (5 minority-spin 
bands vs. 2 majority-spin bands). 

The change in the electronic structure of the one-dimensional NiO chain on the GGA 
level with respect to the B3LYP results can be explained by the insufficient cancellation 
of the self-interaction by the GGA exchange functional, which causes the occupied Ni 3d 
orbitals to artificially rise in energy. 

6.3 0-bridge in Ni nanocontacts 

Atomic chains formed in break junctions have a finite length and are suspended between 
electrodes. It is well known that the contact between the atomic chain and the electrode 
tip will have considerable effect on the electronic structure of the chain, especially when 
d-orbitals are involved like is the case here|55|. 

We have thus calculated the electronic structure and transport properties of both a 
single oxygen atom and a 0-Ni-O chain bridging the two tips of a Ni nanocontact as 
shown in the insets of the right panels of Fig. 16.51 and Fig l6.7l In the case of the single 
oxygen atom (Fig. 16. 5p the electron transport is almost 100% spin-polarized around the 
Fermi level for parallel (P) alignment of the magnetizations of the two Ni electrodes. 
Moreover, an orbital eigenchannel analysis [56| reveals that the transport is mainly due 
to two almost perfectly transmitting minority-spin channels composed of Ni 3di and 
O 2pi orbitals, i.e. they correspond to the doubly-degenerate conduction band of the 
metastable half-metallic state in the perfect chain. Thus the half-metallic state which 
was suppressed in the idealized case of the infinite chain emerges in the more realistic 
situation of a short suspended chain. We can understand this phenomenon in terms of 
the orbital blocking mechanism proposed earlier in the context of Ni nanocontacts [551 [55] . 
The highest minority-spin valence band of the insulating state with FM order in the 
infinitely long chain has a strong contribtution from the Ni 3do orbital which is not 
"compatible" with the geometry of the pyramid shaped Ni contacts, so that this band 
is blocked and thus cannot be occupied. Instead, the doubly-degenerate band composed 
of Ni 3di orbitals hybridized with O 2pi orbitals is partially filled resulting in the half- 
metallic state which in the perfect chain is only metastable. Thus the orbital blocking 
by the contacts actually turns the chain into a half-metallic conductor. Consequently, 
the conductance is strongly suppressed in the case of antiparallel (AP) alignment of the 
magnetizations of the Ni electrodes as can be seen from the right panel of Fig. 16.51 and 
the MR becomes very large: MRi w 90% and MR2 w 700%. 
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Figure 6.5: Transmission per spin channel in the case of P (left) and AP (right) alignment 
of the electrode magnetizations for NiO chain consisting of one oxygen atom bridging the 
two Ni tip atoms of the Ni electrodes as shown on the inset in the right panel. The 
separation of the two Ni tip atoms is 3. 6 A. 



Fig. I6.6f b) shows the transmission per spin-channel for the Ni-O-Ni nanobridge cal- 
culated with GGA. Surprisingly, the transmission does not look very different from the 
transmission calculated with B3LYP for that case. Near the Fermi level the transmission 
is strongly spin-polarizcd: The transmission of the majority-spin channel is strongly sup- 
pressed while in the minority-spin channel essentially two perfectly transmitting channels 
contribute to the conductance. 

This seems to be at odds with the band-structure calculated for the infinite one- 
dimensional chain, which suggests that there should be 5 minority- and 2 majority- 
spin channels contributing to the overall conductance. However, the geometry of the Ni 
nanocontact blocks the transmission of just these channels which due to the unphyscial 
self-interaction of GGA have been raised to the Fermi level. Indeed, an orbital eigen- 
channel analysis [5^ of the transmission reveals that only the doubly-degenerate band 
of type (cZi) contributes to the conductance of the minority-spin channel just as in the 
case of the B3LYP functional. The electrons in the flat minority-spin band of type (^2) 
which actually crosses the Fermi level in the ideal case of the infinite chain are easily 
scattered as they present strongly localized electrons, and thus do not contribute to the 
overall conductance. The other minority-spin channel composed of Ni 3(^3^2 orbitals 
does not contribute either to the conductance since the symmetry of the orbital is not 
compatible with the geometry of the two Ni electrodes - a mechanism to which we have 
referred to in previous work as orbital blocking |55j . The small but finite conductance 
in the majority-spin channel relates to the doubly-degenerate majority-spin band of type 
(di) of the infinite NiO chain raised to the Fermi level due to the self- interaction error. 
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Figure 6.6: (a) GGA band structure of FM groundstate for a lattice spacing of 3.6A. 
Solid lines indicate majority-spin bands and dashed lines indicate minority-spin bands, 
(b) GGA Transmission per for NiO chain consisting of one oxygen atom bridging the two 
Ni tip atoms of the Ni electrodes as shown on the inset in the right panel of Fig. 16.51 The 
separation of the two Ni tip atoms is 3.6A. 
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Figure 6.7: Transmission per spin channel for suspended chain consisting of a O-Ni-0 
bridge and the Ni tip atoms shown in inset of right panel in the case of FM order (left) 
and of AF order (right). For AF order the magnetization of the Ni atom in the center 
of the chain is reversed with respect to the two Ni tip atoms. Distance between a Ni tip 
atom and the center atom is 3.6A. 

The (di) band only crosses the Fermi energy near the upper band edge where the band 
becomes flat. Thus the the electrons in this channel are quite susceptible to scattering 
near the Fermi level resulting in a low transmission. 

Varying the distance d between the Ni tip atoms leads to similiar results as those shown 
in Fig. 16.51 For the P case the current through the chain is almost 100% spin-polarized 
with two open minority channels composed of Ni 3di and O 2pi orbitals, while for the AP 
case it is strongly suppressed, resulting in very high MR values between 80% and 90% for 
MRi for d between 3.0A and 5.0A. Geometry relaxations for different values of the tip-tip 
distance show that for small distances the oxygen atom goes into a zigzag position. The 
bonding angle decreases with increasing distance until it becomes zero at 3.6A. Finally, 
the chain breaks for d > 4.8A. Thus the scattering is strong for small distances d < 3.6A 
when the bonding angle is appreciable and for large stretching, d > 4.2A, resulting in a 
considerable reduction in the conduction of the two minority channels. 

6.4 O-Ni-0 bridge in Ni nanocontacts 

In longer suspended chains the insulating state with FM order starts to emerge inside the 
chain and a away from the contacts. As a result the conductance is reduced as can be 
seen already in the case of the O-Ni-0 bridge (Fig. 16. 7p . The minority-spin conduction is 
reduced considerably (~ 30%) compared to the case of the single oxygen bridge. On the 
other hand the conduction of the majority-spin channel becomes practically zero (< 0.2%). 
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This can be understood by the fact that the residual majority-spin channel conductance of 
the single oxygen bridge is due to direct hopping of Ni s electrons between the electrodes 
and therefore vanishes when the distance between the electrodes is large. The finite 
minority-spin conductance opens up the possibility to an interesting phenomenon: When 
the middle Ni atom reverses its spin, the conductance drops to nearly zero (see right 
panel of Fig. 16. 7p since the AF chain is insulating. In other words, this system behaves 
as a single atom spin valve which presents an extremely large MR even higher than that 
reported above for the single oxygen bridge: MRi 99% and MR2 w 10,000%. Apart 
from controlling the magnetization direction of the central atom by a magnetic field, Fig. 
I6.4f a) suggests that a mechanical control of the spin valve (by stretching the chain) would 
also be possible. 

6.5 Conlcusions 

In conclusion, we have shown that ideal one-dimensional infinite NiO chains could be 
either insulating or halfmctallic depending on the functional used in the calculation. 
Following the results obtained with the B3LYP hybrid functional one-dimensional NiO- 
chains are always insulating for all reasonable values of the lattice spacing but change 
from AF to FM order when stretched slightly out of equilibrium in contrast to bulk NiO 
which has always AF order. A halfmetallic conducting state which is normally metastable 
only becomes the ground state when the chain is stretched to unreasonable large values 
of the lattice spacing. The GGA functional on the other hand predicts NiO chains to 
be ferromagnetic halfmetallic conductors. Although it seems clear that the GGA results 
are affected by the insufficient cancellation of self-interaction error it is a priori not clear 
whether B3LYP predicts the correct ground state as it tends to localize electrons and 
thus tends to predict insulating behaviour. Thus the question of the electronic structure 
of infinite NiO chains remains open at the moment, and more sophisticated many-body 
techniques like GW or DMFT are probably required to answer this question. 

However, in the more realistic case of short NiO chains suspended between Ni nanocon- 
tacts both B3LYP and GGA predict the chains to become strongly spin-polarized con- 
ductors which can be related to the corresponding halfmetallic states of the infinite NiO 
chain, i.e. the metastable halfmetallic state in the case of the B3LYP functional and the 
halfmetallic ground state in the case of the GGA functional. The emergence of almost 
perfect half-metallicity in suspended chains leads to a strong suppression of the current 
for AP alignment of the electrodes resulting in very large MR values of MRi w 90% and 
MR2 ~ 700%, respectively. This could perhaps explain to some extend the very large MR 
values in Ni nanocontacts obtained in some experiments [28|, 135] where oxygen is likely 
to be present. Finally, the 0-Ni-O bridge suspended between Ni electrodes operates as 
a single atom spin valve where the currentflow is controlled by the magnetization of a 
single atom. 
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Fabrication of metallic nanocontacts permits to probe the electronic and mechanical prop- 
erties of conventional metals with unconventional atomic coordination |131] . Electron 
transport in these systems depends on the tiny fraction of atoms in the sample form- 
ing the atom-sized neck which have a reduced coordination and are responsible for the 
two-terminal resistance. Transport experiments can thereby probe the atomic and re- 
lated electronic structure of these atoms and provide information about a fundamental 
question: How bulk properties evolve when the system reaches atomic sizes and atoms 
with full bulk coordination are no longer majority. A bulk property that is susceptible 
to change is magnetism. Bulk Pt, for instance, is a paramagnetic metal but a transition 
to a ferromagnetic state could be expected upon reduction of the atomic coordination 
with the concomitant increase of the density of states (DOS) at the Fermi energy beyond 
the Stoner limit. Density functional calculations [144|, 1145] I146j for one-dimensional infi- 
nite Pt chains support this hypothetical scenario, resulting in a ferromagnetic transition 
above a critical lattice spacing which, depending on the computational approach, can be 
below the equilibrium lattice constant [145] . The formation of local moments in real Pt 
nanocontacts would not be totally unexpected. 

Formation of and electronic transport in finite Pt chains have been extensively studied 
experimentally [1471 11481 11491 11501 1151| . Based upon the appearance of a peak at G = 
0.5 X 2e^//i = 0.5Go in the conductance histogram Rodrigucs et al. suggested that Pt 
and Pd nanocontacts could be spin polarized [1 5 Ij . The origin of this peak has been later 
attributed to adsorbates[35 so that magnetism in Pt and Pd nanocontacts has not been 
confirmed experimentally yet. Previous theory work has addressed the formation of local 
moments in Pd nanocontacts and in Co , Pd and Rh short chains sandwiched between 
Cu planes [135] . To the best of our knowledge theory work on Pt nanocontacts |154l 11551 
I156| I157j has overlooked the possibility of local magnetic order so far. In this chapter 
we perform density functional calculations of both the electronic structure and transport. 
We find that local magnetic order can develop spontaneously in Pt nanocontacts. Local 
magnetic moments as high as 1.2 hi low-coordination atoms are found. Interestingly, 
while transport is definitely spin polarized, the calculated total conductance of magnetic 
and non-magnetic Pt nanocontacts is very similar and in agreement with experimental 
data, explaining why magnetism has been unnoticed so far. 

The electronic structure of various low dimensional structures of Pt which mimic ac- 
tual nanocontacts are calculated in the density functional approximation, using either 
CRYSTAL03[SH] and our ah initio transport package ALACANT that interfaces GAUS- 
SIAN03 to implement the NEGF as explained in based on DFT electronic structure 
calculations Ch. [3l We use scalar relativistic (SR) pseudopotentials for the 60 inner elec- 
trons of the Pt atom and the remaining 18 electrons are treated using generalized gradient 
approximation (GGA) density functionals. The basis set used for all the calculations has 
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been optimized to describe bulk Pt as well as Pt surfaces [158]. Other basis sets such as 
LANL2DZ or SDD[57] have occasionally been employed for comparison. The main results 
do not depend on the choice of basis set. 

7.1 Electronic and magnetic structure of atomic Pt chains 

We first consider a perfect onc-dimcnsional mono-strand Pt chain. Such an idealized 
system serves as a standard starting point to understand lower symmetry geometries. It 
also permits to test whether our LAO pseudopotential methodology reproduces the results 
obtained with SR all electron plane- wave calculations reported by Delin et aL |145| . In 
Fig- IT.lf a) we show the energy per atom as a function of the lattice constant a both for 
the paramagnetic (PM) and the ferromagnetic (FM) chain. They both have a minimum 
at a = 2.4A. The FM chain develops a non-negligible magnetic moment when the lattice 
constant goes beyond a ~ 2.6A. This configuration is clearly lower in energy above that 
distance. The energy difference between the FM and the PM configurations is 16 meV per 
atom for a = 2.7A and 33 mcV per atom for a = 2. 8 A. The magnetic moment per atom 
reaches a saturation value of 1.2^13 ■ The equilibrium distance, critical spacing, asymptotic 
magnetic moment and shape of the phase boundary obtained by us are similar to those 
obtained by Delin et aL |145| using a SR all-electron plane- wave calculation. Our results 
and those of Delin et al. underestimate the onset of the magnetic transition compared to 
calculations including spin-orbit coupling [TJSl [TJS] that predict that a magnetic moment 
forms already below the equilibrium distance. 

The magnetic transition in the phase diagram [Fig. 1(b)] is compatible with the Stoner 
criterion for ferromagnetic instability. As the chain is stretched, the atom-atom coupling 
becomes weaker, the bands narrow down and so does the DOS {T>{e)). Since the integrated 
T>{e) must be equal to the number of electrons per atom, narrowing of the DOS implies an 
increase of the 2?(e) [see Fig. I7.1f d)] and, therefore, an increase of the spin susceptibility, 
which is proportional to 'D{ep). In Fig. I7.1f c) we show how the 'D^ep) of the PM chain 
increases as a function of the lattice constant. The remarkable feature of Pt chains is that 
the Stoner instability occurs close to the equilibrium lattice spacing. 

The electronic structure of the ideal Pt chain sheds some light on the electronic structure 
of the nanocontact. In Fig. 17.21 we show the energy bands for the ideal Pt chain both in 
the FM (left panel) and PM (right panel) configurations, for a lattice spacing a = 2.8A. 
We notice that the spectrum at = has four resolved energy levels per spin. These 
correspond to the 6s level, and the 5d levels which, because of the axial potential created 
by the neighboring atoms, split into 2 doublets E2 and one singlet Ai. The Ei and 
E2 are linear combinations of orbitals with Lz ^ ±1 and = ±2, respectively, whereas 
the Al singlet is a = orbital that hybridizes with the lower energy 6s orbital. The 
largest contribution to the DOS, and therefore to the magnetic instability, comes from 
the Ai-like band at the edge of the Brillouin zone. However, the prominent role played 
by these bands in the magnetic behavior of Pt chains is in stark contrast with their role 
on the transport properties of Pt nanocontacts (see below and see also related work on 
Ni nanocontacts [55]). Four spin-degenerate bands cross the Fermi energy in the PM case 
whereas 7 spin-split bands do it in the FM chain. In the FM chains the number of spin 
minority channels is 6, and the number of spin majority channels is 1. Although spin- 
orbit interaction modifies significantly the bands [145| , the number of bands at the Fermi 
energy is pretty similar in both cases. Therefore, one can anticipate that the number of 
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Figure 7.1: (Color online), (a) Energy per atom for a perfect monostrand Pt chain, (b) 
Magnetic moment per atom as a function of lattice spacing a. (c) DOS of the paramagnetic 
chain at the Fermi energy as a function of a. (d) D.O.S. as a function of energy for 
a = 2.4A (dashed) and a = 2.8A (solid) 
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Figure 7.2: (Color online). Energy 
bands for ideal Pt chain with a = 
2.8A. Left: ferromagnetic phase. 
Right: paramagnetic case. 



k k 

open channels in the magnetic and non-magnetic Pt nanocontacts studied below should 
be roughly the same and thereby the conductance should be similar, although the spin 
polarization might well be large in the former case. The conductance of the ideal FM 
chain is 3.5Go, very far from the value of 0.5Gq that allegedly signals the emergence of 
magnetism and also far away from half of the conductance of the PM chain, so it is very 
unlikely that the celebrated half quantum can be attributed exclusively to magnetism. 

Real Pt chains are typically less than five atoms long and are connected to bulk elec- 
trodes. Although not surprising, we have verified that magnetism survives in isolated 
short chains with Na = 3, 4 and 5 Pt atoms. The equilibrium distance is 2.4A for all 
Na = 3, Na = 4, and Na ~ 5. Interestingly, the short chains are always magnetic in the 
Na = 3 and Na = 4 cases and show a non-magnetic to magnetic crossover at a = 2. 6 A in 
the Na = 5 case, already similar to the ideal infinite chain. The total magnetic moment of 
all the Na = 3 chains with a < 3.0A is 4/iB- The outer atoms have a magnetic moment of 
1.36^B and the central atom with larger coordination has a smaller magnetic moment of 
1.29/iB- In the case of A'^ = 4 the total magnetic moment is 6fiB, and their distribution 
is similar to the Na = 3 case. 

7.2 Transport through suspended Pt chains 

The calculations above show that magnetism is present both in finite- and infinite-sized 
Pt systems with small atomic coordination. It remains to be seen that this holds true in 
nanocontacts where none or only few atoms have a small coordination (like in the case of 
formation of short chains), but these are strongly coupled to the bulk. In order to verify 
whether or not this is the case, we have calculated both electronic structure and transport 
for a model Pt nanocontact. It is formed by two opposite pyramids grown in the (001) 
crystallographic orientation of bulk Pt and joined by one atom which presents the lowest 
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Figure 7.3: Conductance per spin channel for a nanocontact with a 3-atom Pt chain 
(see inset) for the magnetic solution (a) and the non-magnetic one (b). The atom-atom 
distance in the chain is 2.8A. 



possible coordination [see inset in Fig. I7.3f a)]. Relaxation of the 11 inner atoms of the 
cluster has been performed starting from an equilibrium situation as a function of the 
distance d of the outer planes. Zig-zag configurations appear in the chain for small values 
of d (not shown) [T27| until the three-atom chain straightens up (see inset in Fig. I7.3r a)] 
followed by a plastic deformation (not shown). This deformation can be in the form of a 
rupture or a precursor of the addition of a new atom to the chain which comes from one 
of the two 4-atom bases [T44] . 

We now compute the transmission before the plastic deformation occurs (left panel 
in Fig. 17. 3p . where the atom-atom distance in the short chain is 2.82 A. The value 
of the corresponding atomic-plane-averaged magnetic moments are also shown in the 
inset. As expected, it decreases for atoms in the bulk as the coordination reaches the 
bulk value. Some atomic realizations in the stretching process (like zig-zag ones) result 
in nanocontacts with smaller Pt-Pt distance and no magnetism, in agreement with the 
infinite chain phase diagram in Fig. 1. In contrast to the infinite chain, there are only 
three channels contributing to the total conductance for minority electrons and there are 
more than one (three) for majority ones. For the majority electrons these are a perfectly 
transmitting s-type channel and two partially open (T = 0.4) pd-type channels {one Pxdxz- 
and one Pydyz- hybridized). The three minority channels have the same character as the 
majority channels, except that here the s-type channel does not transmit perfectly while 
the transmission of two pd-type channels is enhanced so that all three minority channels 
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have a transmission around 0.7. The other two remaining pd-hke channels are responsible 
for the sharp resonances that appear around the Fermi level. The total conductance 
of the nanocontact in the FM case thus turns out to be around 4e^//i = 2Go which is 
only slightly larger than the average experimental value corresponding to the last plateau 
(1.75Go), but, interestingly, barely differs from the value obtained when the possibility 
of magnetic order is ignored [around 2.36*0, see right panel in Fig. 17. 3j . As in the case 
of the ideal chain, the FM conductance is not half of the PM conductance nor half of 
Gq. To conclude this discussion we notice that although transport is only weakly spin 
polarized, magnetism brings the pd-like resonances up to the Fermi level compared to the 
non-magnetic case. These resonances may well give features in the low bias conductance 
not present if Pt were not magnetic. 

7.3 Conclusions and Discussion 

The main conclusion of this part is that density functional calculations predict that 
nanochains formed in Pt nanocontacts can be stretched as to become magnetic. The 
magnetic moment is localized mainly in the atoms with small coordination and does not 
modify appreciably the total conductance, although the transmission is moderately spin 
polarized. How robust are these results? It is well known that both local and gradient- 
corrected density functionals present some degree of electronic self-interaction, in contrast 
with the Hartree-Fock approximation. Self-interaction is larger for localized electrons and 
shifts the d bands upwards in energy, as shown in the case of Co, Ni and Pd one dimen- 
sional chains p"59] . A number of schemes to avoid this problem, like LDA-I-U and Self- 
Interaction Correction functionals have been proposed. The method of choice between 
chemists is hybrid functionals (51] in which local and Hartree-Fock exchange (HFX) are 
combined and the self-interaction is reduced. We have calculated the magnetic phase di- 
agram of the one dimensional Pt chain using the hybrid B3LYP functional [SI] and found, 
somewhat expectedly, that magnetism is enhanced and that B3LYP infinite Pt chains are 
ferromagnetic down to the equilibrium distance (a = 2.4A). Both non-local exchange and 
spin-orbit coupling [145] enhance the stability of magnetism in Pt nanocontacts. This 
and previous results on Ni nanocontacts [5S] lead us to believe that self-interaction is an 
issue in the electronic structure and transport properties of transition metal nanocontacts 
and further work is necessary along these lines [160] . 

The mean field picture of the electronic structure describes a static magnetic moment 
without preferred spatial direction. In reality the nanomagnet formed in the break junc- 
tion is exchanged coupled dynamically to the Fermi sea of the conduction electrons of 
the electrodes. In the case of a spin 5 = 1/2 this can result in the formation of a Kondo 
singlet that would yield an anomaly in the zero bias conductance. For larger spins, the 
conduction electron sea cannot screen the spin completely so that the magnetic moment 
survives. The magnetic moment of the nanocontact in Fig. 3 is 5* w 6, comparable to 
that of single molecule magnets |161] , making the formation of a Kondo singlet unlikely. 

In the absence of spin-orbit interactions and external magnetic field a electronic config- 
uration with total spin S has 25* -I- 1 degenerate configurations corresponding to the spin 
pointing along different directions. However, spin-orbit interaction is strong in Pt and pro- 
duces spin anisotropy, favoring orientation along the transport direction axis in the case of 
one dimensional chains [25] ■ Thermal fluctuations of the magnetic moment between these 
two configurations arc quenched for temperatures smaller than the anisotropy barrier. 
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Departures from the easy axis orientation will be damped via electron-hole pair creation 
across the Fermi energy that would also result in small bias features in transport [162] . 
The application of a sufficiently strong magnetic field in the direction perpendicular to 
the easy axis moves the local magnetic moments away from their easy axis. This is known 
to change the number of open channels at the Fermi energy in both the case of Ni ideal 
chains [163] and in the case of ferromagnetic semiconductor tunnel junctions [164j . This 
effect, or maybe even larger, can be expected in Pt nanocontacts and could be used to 
detect the nanomagnetism experimentally. 
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8.1 Overview of developed computational tools 

In this thesis, I have studied spin transport through nanocontacts and nanowires, both 
with calculations of simplified models and with full ab initio calculations based on density 
functional theory (DFT). To this end various computational tools have been developed 
either from scratch or by extension of already existing code: 

i) SpinTrans: This program for computing spin transport in simple models with 
Coulomb interaction (e.g. Hubbard model) in the non-collinear unrestricted Hartree- 
Fock approximation (NC-UHF), see also ApplEl has been developed from scratch in the 
first stage of the thesis. This program was applied for the self-consistent calculation of 
toy models of magnetic nanostructures to study the effect of non-collinear magnetization 
profiles on the transport presented in Ch. IH 

ii) ALACANT: The ALACANT (ALicante Ab-initio Computation Applied to Nano 
Transport) project was started in the year 2000 by professor J. J. Palacios in the Applied 
Physics department of the University of Alicante. It implements the NEGF formalism 
in connection with ab initio electronic structure calculations based on DFT. as explained 
in ChlJl As a part of this thesis, the ALACANT package has been further developed in 
various aspects: 

• The spin-unrestricted transport formalism was implemented into the package 

• The entire program code which had been written in FORTRAN?? has been ported 
to FORTRAN90 in order to make use of modern programming techniques like the 
use of modules, dynamic memory allocation etc. 

• A module for calculating self-energies for one-dimensional electrodes described by 
Hamiltonian and overlaps matrices which can be taken from ab-initio calculations 
has been implemented. 

• An interface to the CRYSTAL ab initio program for crystalline systems has been 
developed. 

8.2 Summary of results 

With the developed computational tools various calculations of the electronic and mag- 
netic structure and the transport properties of nanocontacts and nanowires have been 
performed. First, we have studied spin transport with simple models, namely a one- 
dimensional Hubbard chain in the ferromagnetic phase, in order to gain an understand- 
ing of the fundamental mechanisms of nanoscale spin transport. Neglecting the geometric 
aspects as well as the much more complex Hamiltonian and Hilbcrt space of real nanocon- 
tacts allowed to concentrate on the pure spin aspects of electron scattering by the magnetic 
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structure. We find that domain walls form self-consistently under the appropriate mag- 
netic boundary conditions. Furthermore, we find that the longer the atomic chain the 
smoother the domain wall becomes, and the scattering by the domain wall is reduced 
since the spin-mixing transforms the spin of an incoming electron adiabatically into the 
opposite spin. Vice versa, for short necked nanocontacts (i.e. no chain formation) like is 
typical for Ni, the domain walls will be rather sharp, and thus the spin scattering by a 
domain wall formed in the neck of the nanocontact becomes maximal. 

In order to investigate the question of the possibly huge BMR values in Ni nanocon- 
tacts we have performed ab initio calculations of the electronic structure and transport 
properties of Ni nanocontacts on the level of DFT. By comparing solutions with and 
without formation of a domain wall in the atomic neck of the nanocontact, we find that 
BMR is certainly not large in pure Ni nanocontacts, but rather moderate. This is due to 
the fact that a grand part of the a priori available spin-polarized channels which could 
in principle give rise to a large BMR are blocked by the geometry of a real nanocontact. 
On the other hand the spin-unpolarized s-type channel which does not give rise to any 
BMR is not affected by the geometry. Thus even in the DW configuration the conduc- 
tance is appreciable 2e^/h compared to the < 4e^ for the ferromagnetic solution). 
Another important conclusion that we can draw from our calculations is that disorder of 
the atoms in the nanocontact strongly reduces the spin-polarization of the current since 
it has a stronger effect on the transmission of the spin-polarized d-typc channels than on 
the unpolarized s-typc channels. Moderate BMR values for ferromagnetic nanocontacts 
has been confirmed recently by experiments with very clean samples under controlled con- 
ditions (ultra-high vacuum conditions, exclusion of any magnetostriction effects) |38l 137] . 
and also by a number of theoretical papers P^5] . 

Next, we have studied the electronic structure and transport properties of one-di- 
mensional NiO chains, both ideal infinite ones, and short ones suspended between Ni 
nanocontacts. Indeed, it has been observed experimentally, that Ag which like Ni does 
not form atomic chains in a nanocontact. actually docs form chains when oxygen is present 
[132j . and it has been argued that the incorporation of oxygen atoms into the a chain 
actually stabilizes the chain |165j . In fact we find that NiO chains bridging the tip atoms 
of a Ni nanocontact are stable, and so NiO chains could in principle form when nanocon- 
tacts are fabricated in an oxygen atmosphere. Most importantly, we have found that 
short NiO chains can actually become almost perfect half-metallic conductors when sus- 
pended between Ni nanocontacts, i.e. the current becomes almost 100% spin-polarized. 
In fact, already a single oxygen atom bridging the tips of a Ni nanocontact can have such 
an effect. Consequently, this leads to a huge BMR of the order of the GMR effect, and 
could perhaps explain to some extent the large MR values obtained in some experiments 
[281 135j . Another very important conclusion we can draw is that already a single atom 
can completely change the conductance behavior of an atomic-scale device. 

The emergence of magnetism at the nanoscale in materials that are otherwise param- 
agnetic in bulk is a fascinating topic of Nano science, and has been observed for example 
in Au nanoclusters |166| . and predicted for atomic chains of Pt and Pd with ab initio 
electronic structure calculations [SJl [55] . In order to see whether conduction measure- 
ments of Pt nanocontacts could possibly reveal the magnetism of an atomic chain formed 
between the tip atoms of the nanocontact, we have performed ab initio electronic struc- 
ture and transport calculations of atomic Pt chains suspended between the tips of a Pt 
nanocontact. We find that although suspended Pt chains are also magnetic, the conduc- 
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tance is very similar to the case when there would be no magnetism in the chain. Thus it 
seems not possible to decide experimentally whether Pt chains become magnetic or not 
by simple conductance measurements of Pt nanocontacts. 
Furthermore, we can draw the following general conclusions: 

i) The DFT based transport approach on the LDA level at small bias voltages seems to 
work reasonably well for metallic nanocontacts and nanowires, even when d-electrons are 
involved in the conduction. This is in agreement with theoretical studies which show that 
DFT calculations on the LDA or GGA level give a reasonable description of the electronic 
structure of bulk metals and metallic surfaces |167l 1115] , and corrections by methods like 
the GW approximation or DMFT which take into account electron correlations are often 
relatively small in the vicinity of the Fermi level |168|, 1169] . However, when going away 
from the purely metallic systems things become more complicated. NiO for example is 
a strongly correlated material, and both LDA and GGA fail in describing its electronic 
structure correctly. Also the problems of LDA/GGA of describing the electronic structure 
of molecules are notorious. Thus it must be doubted whether LDA or GGA can give a 
good description of nanoscale NiO e.g. oxidized Ni nanocontacts or one-dimensional NiO 
chains, or molecular conductors. On the other hand the semi-empirical DFT methods (e.g. 
hybrid functionals or LDA-I-U) often describe quite well the electronic structure of some 
bulk oxides with strongly correlated electron systems, and especially hybrid functionals 
are extremely successful in describing molecules. However, these methods require an 
adjustment of the empirical parameters by fitting with experimental results, which limits 
their predictive power. Also hybrid functionals often do not yield a satisfying description 
of metallic systems. Thus the electronic structure of a molecule in a molecular electronics 
device probably needs to be described by a hybrid functional while the metal electrodes 
for contacting the molecule are not well described by hybrid functionals but needs to be 
described on the LDA or GGA level. 

ii) Magnetic nanocontacts and nanowires are quite sensitive to the atomic structure, 
i.e. distortions and defects often have a strong effect on the transport properties. This 
can be understood by the fact that the d-orbitals which are in fact responsible for the 
magnetism are very sensitive with respect to the geometry and thus d-electrons are easily 
scattered. On the other hand, the conduction channels composed of d-orbitals are the 
only channels that are actually spin-polarized and thus are the ones responsible for the 
magneto-resistive behavior. So exactly the channels that are interesting in the context of 
spintronics applications are not stable against variation of their atomic structures. Thus it 
is arguable whether magnetic atomic-size nanocontacts could one day serve as ingredients 
for nanoscale electronics. But rather they allow us to gain a fundamental understanding 
of conduction and the interplay with other physical processes at the atomic scale, and to 
develop an appropriate methodology for the description of nanoscale conductors. This 
might well proof essential for the development of new electronic devices in the near future. 

8.3 Outlook 

For the results presented in this thesis the spin-orbit coupling of the electrons has been 
neglected. This can be justified by the fact that the spin-orbit coupling is relatively small 
for not too heavy elements like Ni. For heavier elements like Pt on the other hand the 
spin-orbit coupling becomes comparable to the other energies and thus probably affects 
the electronic and magnetic structure of the material: As pointed out at the end of Ch. [7] 
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8 Summary and Outlook 

the spin-orbit interaction introduces spin-anisotropy favoring the magnetization of the Pt 
chain to be oriented along the transport axis. Furthermore the spin-orbit couphng favors 
the onset of magnetism in atomic Pt chains as has been shown by ab initio calculations 
including the spin-orbit coupling [145^ I146j . 

By applying a sufficiently strong magnetic field, the orientation of the magnetization 
of the chain atoms can be changed from the one favored by the spin-orbit coupling called 
the easy axis. This however affects the electronic structure of the atomic chain due to 
the coupling of the orbital degree of freedom to the spin degree of freedom and possibly 
alters the transport properties of the system resulting in anisotropic magneto-resistance 
(AMR) defined as the maximal change in resistance when rotating the magnetic field |163j . 
Despite the relatively small spin-orbit coupling of Ni and Fe in the order of 50meV, recent 
experiments have measured an appreciable AMR of 10-15% in the case of Ni nanocontacts 
[551 [T7D) and huge AMR of up to 75% in the case of Fe nanocontacts [TTT] . 

We have therefore implemented spin-orbit coupling into our ab initio transport program 
ALACANT and obtained first results for one-dimensional Ni chains with a scattering 
center as simplified models of real Ni nanocontacts. The calculated AMR values are 
similar to the experimentally measured AMR values. However, further investigations with 
more realistic geometries are needed in order to compare with the experiments as the spin- 
orbit coupling solely affects the d-electrons which are very sensitive to geometrical effects 
as has been discussed extensively in this thesis. Due to the strong spin-orbit coupling 
of Pt one might also expect that the magnetism of atomic Pt chains should give rise to 
a strong AMR signal which would thus give experimental evidence for the emergence of 
magnetism in atomic Pt chains. Thus it could be worthwhile to calculate the AMR of 
Pt nanocontacts from ab initio quantum transport calculations including the spin-orbit 
coupling. 

An important problem for future nanoscale spintronics applications is the injection 
of spin-polarized currents from magnetic into non-magnetic nanoscopic conductors like 
paramagnetic nanowires, carbon nanotubcs and organic molecules. The possibility of 
spin-injection into a non-magnetic material will very likely depend strongly on the atomic- 
scale properties of the contact or between the magnetic and non-magnetic material. It is 
thus of fundamental importance to study interfaces between magnetic and non-magnetic 
nanoscale conductors. Related to the question of spin injection is the problem of spin re- 
laxation which designates the phenomenon that the spin-polarization of a current through 
a non-magnetic material decreases due to the precession of the electron spin caused by 
the spin-orbit coupling. Studying spin injection and spin relaxation in nanowires, nan- 
otubes and organic molecules with ab initio quantum transport calculations thus defines 
an interesting future line of work. 
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A Representation of operators in 
non-orthogonal basis sets 



The natural definition for the matrix A of a one-body operator A in a non-orthogonal 
basis set (NOBS) {|a)} is simply by its matrix elements: 

A= = . (A.l) 

However, the representation of an operator in a NOBS is not that simple; 

i = ^|a)(S-iAS-i)„^3(/3| (A.2) 

where S = [Saf^) = is the overlap matrix for the basis set. It is easy to see that 

this definition leads results in the matrix elements Aafj defined above. Then the identity 
operator in the NOBS representation is given by 

i = ^|a)(S-i),/3(/3|, (A.3) 

which is also easy to proof. 

Now we define a second matrix 

A:=S-iAS-\ (A.4) 

which is the matrix that appears above in the representation of the operator in a NOBS. 

One should take care when using the representation of an operator in a NOBS. For 
example, the matrix element of an operator between two non-orthogonal orbitals |a), \0) 
can be zero, (q;|A|/?) = 0, but the corresponding matrix element of the matrix A does 
not necessarily vanish due to the multiplication with the inverse of the overlap matrix 
on both sides. Thus there is actually a non-zero contribution of the two orbitals to the 
operator although the corresponding matrix element of the operator is zero 

Orthogonalizing the basis set by the Lowdin orthogonalization scheme [70] , the matrices 

A and A are transformed to the matrix A-*- ~ defined in the new orthogonal 

basis set {|«)} according to: 

A^ = S-1/2AS-1/2 = (A.5) 

Though there are also other orthogonalization schemes, the Lowdin scheme is particularly 
useful in the context of quantum chemistry methods based on atomic orbitals as the 
center of the orthogonalized orbital remains centered on the same atom as the original 
non-orthogonal orbital. 
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A Representation of operators in non-orthogonal basis sets 
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B Partitioning method 



As explained in Ch. [5] we model the transport problem by dividing the system in three 
parts. Two semi-infinite leads (L) and (R) with bulk electronic structure are connected 
to a finite region called device (D). In a local basis set the Hamiltonian and the overlap 
matrix of the system are given by (|2.13p and (|2.14[) . Dividing the F matrix into sub- 
matrices in a similar manner we obtain the following matrix equation: 



zSl-Hl 


zSld 


— Hld 







z Sd 


He 


zSdr 




zSrd 


— ii-RD 


zUr 




Gl{z) Gj^Diz) GLRiz) 
Gdl{z) Goiz) GDfl(z) 
GrUz) Grd^z) Gfl(z) 




(B.l) 



This yields 9 equations for the 9 sub-matrices of the GF G. We can resolve this matrix 

equation columnwise. Multiplying all rows of ES — H with the first column of G yields 
three equations for Gl, Gjjl and Grl which yield: 

Gl{z) = {zSi^-UL-f:D+R{z)r^ (B.2) 

Gdl{z) = gD+R{z){liDL-zSDL)GL{z) (B.3) 

Grl{z) = SR{z){ilRD-zSRD)GDL{z) (B.4) 
Similarly we obtain from multiplication with the second column: 

Gd{z) = {zSD-liD-^L{z)-^R{z))-^ (B.5) 

GLoiz) = gL{z){liLD-zSLD)GDiz) (B.6) 

GRoiz) - sr{z){11rd-zSrd)Gd{z) (B.7) 
And finally from multiplication with the third column, we obtain: 

Gr{z) = {zSr-IIr-^d+l{z))-^ (B.8) 

GdrXz) - gD+L{z){llDR-zSDR)GR{z) (B.9) 

Glr{z) = SL{z){liLD~zSLD)GDR{z) (B.IO) 

We have introduced the Green's functions of the isolated left and right lead and 
and the corresponding self-energies El and T,r: 



sdz) 

^l{z) 
Sr{z) 



{zSl-Hl)-' (B.ll) 

(Udl - zSDL)gL{z) {Hld - zSld) (B.12) 

(zSr-Hr)-^ (B.13) 

{B.DR-zSDR)gR{z)iURD~zSRD) (B.14) 
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B Partitioning method 

Furthermore, we have defined the Green's function of the device plus the left lead only, 
gD+L, of the device plus the right lead only, go+R, and the corresponding self-energies 
S D+ L and S D+R each one representing the coupling of one of the leads to the device and 
the other lead: 



^d+l{z) 

gD+R{z) 



{zSd-^d-^l{z))-^ 
(zSz5-Hi3-Efl(z))-i 

(H_R_D - zSiiD)gD+L{z){^DR - zSdr) 
(HlZ5 - zSLD)'iD+R{z){'iiDL - zSdl) 



(B.15) 
(B.16) 
(B.17) 
(B.18) 
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C Self-energy of a one-dimensional lead 



Here we will derive the Dyson equation (|2.42p for the calculation of the self-energy of the 
semi-infinite right lead. The derivation of the Dyson equation for the left lead ()2.4ip goes 
in a completely analogous way. 

The Hamiltonian matrix Hn of the (isolated) semi-infinite right electrode is defined in 
eq. (PTTB)) as: 



H 



R — 



H, 



H Hn 



Hi 



Hi 
Ho 



Hi 



Vo 



(C.l) 



and the overlap matrix is given in eq. (j2.18p as: 



Vo 



Si 

50 Si 

51 So Si 



(C.2) 



To obtain the self-energy of the lead we have to calculate the GF of the lead from its 
defining equation: 

{zSr - Hfl)iH(z) = 1. (C.3) 

In the same way as the Hamiltonian and the overlap matrix we subdivide the GF matrix 
gfl, into sub-matrices corresponding to the unit cells of the lead. Now the above equation 
for the right lead's GF reads: 



/zSo-Ho 2S1-H1 



zS 



H\ zSq - Ho zSi - Hi 









1 



■ / 

(C.4) 



As explained in Sec l2.2l it suffices to calculate the "surface" GF, i.e. gi.i. From multi- 
plication of the 1st, the 2nd and so on until the n-th line of {zSn — Hr) with the 1st 
column of g_R(^) we get the following chain of equations: 

(zSo-Ho)gi,i(2) + (^Si-Hi)g2,i(^) = 1 (C.5) 
(zSl-Hl)gia(z) + (^So-Ho)g2,i(;^) + (zSi-Hi)g3,iW = (C.6) 



(zSl-Hl)g„_i,i(z) + (zSo-Ho)g„,i(z) + (zSi-Hi)g„+i,i(z) = (C.7) 
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C Self-energy of a one-dimensional lead 



For n> 1 the equations for determining gn,i{z) all have the same structure: 

(0So-Ho)i„4(z) = (Hi -zSl)g„_i,i(z) + (Hi -zSi)g„+i,i(z). (C.8) 

We define a transfer matrix for n > 1 by: 

T„_i,„(z)g„_ia(z) = in,l(2). (C.9) 

The transfer matrix thus transfers information from site n — 1 to site n of the lead, i.e. 
from the left to the right. Multiplying Eq. (|C.8p by (g„_i^i)~^ we obtain: 

(zSo - Ho) T„-i,„(z) = (Hi - zS\) + (Hi - zSi) T„,„+i(z) T„_i,„(z) (C.IO) 

Reordering we obtain the following iterative equation for the transfer matrices: 

T„_i,„(z) = (zSo - Ho - (Hi - zSi) T„,„+i(z))"i - zS{) (C.ll) 

Since the electrode is semi-infinite it looks the same from each unit cell when looking to 
the right. Thus a given g„_i^i, results always in the same g„^i independent of n. Thus 
the transfer matrix must be independent of n: T„_i^„(z) = T(z), and Eq. (|C.11[) allows 
to determine the T(z) self-consistently. 

We define the self-energy as := (Hi — zSi) T(z), and obtain the Dyson equation 
for the self-energy: 

= (Hi - zSi) (zSo - Ho - S(z))-' - zS\). (C.12) 

We will now see that this self-energy is indeed identical to the one defined for the right 
lead in eq. ()2.38p . i.e. S(2:) = S,.(£'). By plugging in the definition of the transfer matrix, 
eq. (|C.9p . into eq. (|C.5|) for determining the surface GF, gi_i we find: 

(zSo-Ho)gia(z) + (zSi-Hi)T(z)gi,i(z) = 1 

^ (zSo-Ho-t-S(z))gi,i(z) 1, (C.13) 

where in the last step we have made use of the definition of the self-energy. Thus we 
obtain for the surface GF of the right lead: 

gi,i(z) = (zSo-Ho + S(2))-i. (C.14) 

And vice- versa the self-energy can be expressed in terms of the surface GF: 

S(z) = (Hi - zSi) gi,i(z) (Hi - zS\). (C.15) 

This proofs that the self-energy S(z) defined above in terms of the transfer matrix is 
identical to the self-energy Sr(z) defined earlier in Sec l2.2l so that the self-energy Tir{z) 
can be calculated iteratively by the Dyson equation (|2.42p . 

The proof for the left lead runs completely analogously. The surface GF of the left lead 
is now: 

= (zSo-Ho + S,(z))-i. (C.16) 
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D Bethe lattices 



In this appendix we discuss how self-energies for Bethe lattices (BL) used to describe the 
leads are calculated. A BL is generated by connecting a site with N nearest-neighbors 
in directions that could be those of a particular crystalline lattice. The new N sites are 
each one connected to TV — 1 different sites and so on and so forth. The generated lattice 
has the actual local topology (number of neighbors and crystal directions) but has no 
rings, and thus does not describe the long range order characteristic of real crystals. Let 
n be a generic site connected to one preceding neighbor n—1 and N —1 neighbors of the 
following shell {n + i with i — 1, .., N— 1). Dyson's equation for an arbitrary non-diagonal 
Green's function is 

— Ho)G„^fc = V„^„_iG„_i^fc + ^ 'Vn.iGiM (D.l) 

i=l,...,JV-l 

where k is an arbitrary site, E the energy, and V^.j is a matrix that incorporates in- 
teractions between orbitals at sites i and j (bold capital characters are used to denote 
matrices). Hq is a diagonal matrix containing the orbital levels and I is the identity 
matrix. Then, we define a transfer matrix as 

Ti_i_iGi_ij- = Gij (D.2) 

Multiplying Eq. (|D.1[) by the inverse of G„_i.„ we obtain, 

(i?I — Ho)T„_i^„ = V„.„_i + I ^ ^n,iTn,i I T„_i^„ (D.3) 

= l,...,Ar-l 

Due to the absence of rings the above equation is valid for any set of lattice sites, and, 
thus, solving the BL is reduced to a calculation of a few transfer matrices. Note that a 
transfer matrix such as that of Eq. (|D.2p could also be defined in a crystalline lattice but, 
in that case it would be useless. 
Eq. (jP.Sp can be solved iteratively, 

-EI - Ho - J2 ^r,,^T„, 
.N-1 



L n— l,n 



Vn,n-1 (D.4) 



If the orbital basis set and the lattice have full symmetry (including inversion symmetry) 
the different transfer matrices can be obtained from just a single one through appropriate 
rotations. However this is not always the case (see below). 

Before proceeding any further we define self-energies that can be (and commonly are) 
used in place of transfer matrices, 

S,,, = V,,,T,j (D.5) 
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D Bethe lattices 



Eq. (jD.4|l is then rewritten as, 



l.n "^Ti— l,n 



SI — Ho ^ S„,i 



i=l,...,JV-l 



(D.6) 



where we have made use of the general property Vn,n-i = v]^_j „. 

As discussed hereafter, in a general case of no symmetry this would be a set of N coupled 
equations (2A'^ if there is no inversion symmetry). Symmetry can be broken due to either 
the spatial atomic arrangement, the orbitals on the atoms that occupy each lattice site, 
or both. When no symmetry exists, the following procedure has to be followed to obtain 
the self-energy in an arbitrary direction. The method is valid for any basis set or lattice. 
Let Ti be the N nearest-neighbor directions of the lattice we are interested in and Vr-^ the 
interatomic interaction matrix in these directions. To make connection with the notation 
used above note that the vector that joins site n — 1 to site n, namely, rn — r„_i would 
necessarily be one of the lattice directions of the set Tj. The self-energies associated to 
each direction have to be obtained from the following set of 2N coupled self-consistent 
equations. 



S,^ = V,^ [EI - Ho - - ^nT 



^r, = Vn [EI - Ho - (St - SrJ]"' V 



t 



(D.7) 
(D.8) 



where i = 1,...,N and fi = —Ti. V^-. is the interatomic interaction in the direction, 
and St and are the sums of the self-energy matrices entering through all the Cayley 
tree branches attached to an atom and their inverses, respectively, i.e.. 



N 



(D.9) 



N 



(D.IO) 



This set of 2N matricidal equations has to be solved iteratively. It is straightforward to 
check that, in cases of full symmetry, it reduces to the single equation. The local density 
of states can be obtained from the diagonal Green matrix. 



EI Ua 

i=l,..,N 



(D.ll) 
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E Non-Collinear Unrestricted 
Hartree-Fock Approximation 



The standard Unrestricted Hartree-Fock theory (UHF) only allows for spin-polarized 
ground states, i.e. the UHF Hamiltonian commutes with one component of the total 
spin operator, e.g. Sz, so that the ground state j^'o) is an eigenstate to Sz- For atoms 
and molecules this is a good property of the UHF Hamiltonian. However, in transport 
problems, this symmetry might be broken if the magnetizations of the two macroscopic 
leads are aligned non-coUinear or anti-parallel. 

Therefore, we will derive a non-coUinear formulation of unrestricted Hartree-Fock the- 
ory (NC-UHF). As in UHF we want to minimize the total energy of a single Slater- 
determinant of an A'^-electron system, 

|*o) = |xi,...,Xw) (E.l) 
described by the many-body Hamiltonian 

i i<j 

Only now the molecular orbitals are a linear combination of atomic orbitals with different 
spins, i.e. we explicitly allow for spin-mixing: 

i.cr 

where |0f ) = ® \a). 

We assume that the molecular orbitals and the atomic orbitals are orthonormal: 

{Xa I Xb) = Sab and {(f>i \ (pj) = % (E.4) 

In order to find the molecular orbitals that minimize the total energy with the con- 
strained of orthogonal molecular orbitals we have to minimize the Lagrangian 

>C[{Xa}] = Eo[{Xa}] - J2 iiXa \ Xb) - Sab) , (E.5) 

a, 6 

where Eo[{xa}] = (*o|7t;|*o). 

From the minimization condition SC = and the expression for the total energy of a 
single Slater determinant 

Eo = ^{Xa\hQ{l)\xa) + '^^{Xa,Xb\v{l,2)\xa,Xb) 
a a.b 

- (^-'^"1^(1' 2)1X6, Xa) (E.6) 

a, 6 
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E Non-Collinear Unrestricted Hartree-Fock Approximation 



wc find the following effective single particle Hamiltonian - the Fock operator 

i,j-cr k,l 

Xl^;fe'''(^''''^fc|*(l'2)|0(,?!ij)cJ^^Cj^, (E.7) 

i,j,cri,(T2 k,l 

which depends on the density operator 

a i,j.(7i ,(72 a— 1 

The second term in cq. (jE.7[) is the Hartree-term -O/f which describes the classical 
(direct) Coulomb repulsion of the electrons. It depends only on the (local) charge density 

VH^ j dx j Ayp^{y)v{x, y) E '4^l{x)'ipa{x), (E.9) 

while the third term really depends on the density matrix. This is the Fock-tcrm vp 
describing the exchange energy: 

VY^Y. J'^'^J dyp'''^Hx,yHx,y)^Ui^)'P-2iy)- (E.io) 

Clearly, only the Fock-term can give rise to non-coUincar spin-arrangements. Moreover, 
only if the initial guess contains non-diagonal terms with respect to the spin, i.e. only 
if p'^^"'^{x,y) for cti (J2, the self-consistent Hartree-Fock solution can have a non- 
collinear spin configuration. 
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